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Abstract 

We study a quasiperiodic structure in the time evolution of the Bloch vec- 
tor, whose dynamics is governed by the thermal Jaynes-Cummings model (JCM). 
Putting the two-level atom into a certain pure state and the cavity field into a mixed 
state in thermal equilibrium at initial time, we let the whole system evolve according 
to the JCM Hamiltonian. During this time evolution, motion of the Bloch vector 
seems to be in disorder. Because of the thermal photon distribution, both a norm 
and a direction of the Bloch vector change hard at random. In this paper, taking 
a different viewpoint compared with ones that we have been used to, we investi- 
gate quasiperiodicity of the Bloch vector's trajectories. Introducing the concept of 
the quasiperiodic motion, we can explain the confused behaviour of the system as 
an intermediate state between periodic and chaotic motions. More specifically, we 
discuss the following three facts: (1) We show that the dynamics of the Bloch vec- 
tor is equivalent to that of a compressible fluid. (2) If we adjust the time interval 
At properly, figures consisting of plotted dots at constant time intervals acquire 
scale invariance under replacement of At by sAt, where s(> 1) is an arbitrary real 
but not transcendental number. (3) We can compute values of the time variable 
t, which let (S^i)] (the absolute value of the z-component of the Bloch vector) be 
very small, with the Diophantine approximation (a rational approximation of an 
irrational number). 



1 Introduction 

The Jaynes-Cummings model (JCM) describes a two-level atom interacting with a single- 
mode radiation field. Because it is a soluble fully quantum mechanical model, the JCM has 
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been attracting many researchers' attention since it was proposed in 1963 [1, 2, 3, 4, 5]. 
Considering the Hamiltonian for a magnetic dipole in a magnetic field, assuming near 
resonance and applying the rotating- wave approximation to it, we obtain the Hamiltonian 
of the JCM. Thus, the JCM is suitable for describing interaction of the radiation with the 
matter, for example, the spontaneous emission given rise to by a single two-level atom in 
the cavity field. 

If we put the single-mode field into the state with sharply defined photon number 
at initial time, the JCM causes the Rabi oscillations in the atomic population inversion. 
Moreover, if we put the field into the coherent state initially, the JCM shows the collapse 
and the revival of the Rabi oscillations. Because the classical theory cannot explain the 
revival of the Rabi oscillations, we can regard it as an evidence of the quantum nature of 
the radiation, that is, the discreteness of photons [6, 7, 8]. 

However, in the case where we prepare the single-mode field in a thermal state initially, 
we hardly predict how the system evolves. Because of the thermal fluctuation of the field, 
the Bloch vector develops in time in a confusing manner. Both its norm and direction 
change hard at random, so that it seems to be in disorder. A number of researchers 
analyze this problem in details and try to find distinct properties, which characterize its 
unpredictable behaviour. 

In Ref. [9], von Foerster computes time evolution of the probability for finding the 
atom in the ground state with preparing the radiation field into various states at ini- 
tial time. Knight and Radmore investigate both the thermal JCM and its corresponding 
semiclassical model [10]. They compare time evolution of the atomic population inver- 
sion of the quantum model with that of the semiclassical model by plotting their graphs 
numerically. In Ref. [11], Knight analyzes the time-dependence of the atomic inversion 
of the JCM with an initial field state in thermal equilibrium. In the work, he states the 
following observation: The thermal fluctuation causes the collapse of the Rabi oscillations 
and their collapse time is very short. In addition, the revival overlaps it, so that inter- 
ference occurs and it generates a very irregular time evolution. Liu and Tombesi discuss 
thermodynamics of the JCM [12]. In their analysis, the grand partition function of both 
the atom and the boson field is considered. 

When we discuss the thermal JCM, we have to handle an intractable infinite series. If 
the single-atom field is resonant with the atom, the nth term of this intractable infinite 
series is given by a trigonometric function of \J~nt for n — 0, 1,2, where t represents 
the variable of the time. Buck and Sukumar emphasize this fact in Ref. [13]. Because the 
infinite series is a superposition of the trigonometric function of ^fnt for n — 0, 1, 2, ... and 
it cannot be a Fourier series, the value of the sum of the series varies in an unpredictable 
manner as the time t progresses. Arroyo-Correa and Sanchez-Mondragon try to discuss 
the thermal JCM and evaluate the atomic population inversion, which is described by 
this intractable infinite series, using a technique of the complex analysis [14]. Chumakov 
et al. examine a new analytical approach to obtain an approximate sum of this series, 
which is reliable for a small initial mean photon number [15]. In Ref. [16], Klimov and 
Chumakov rewrite this intractable infinite series as a sum of two integrals by using the 
Abel-Plana formula. 

Not only the thermal fluctuation but also a continuous measurement lets the JCM 
exhibit chaotic behaviour. Fukuo et al. study time evolution of non rotating-wave ap- 



2 



proximation JCM under a continuous quantum-nondemolition (QND) measurement [17]. 
Counting the number of photons inside the cavity with the QND photodetector, the sys- 
tem evolves in an irreversible way and yields chaos. 

In this paper, we study the time evolution of the Bloch vector governed by the ther- 
mal JCM. Preparing the atom and the single- mode field initially in a certain pure state 
and a mixed state with a Bose-Einstein photon number distribution respectively, we let 
the whole system develop according to the JCM Hamiltonian. Thus, the system never 
suffers from dissipation. As mentioned before, under these circumstances, the behaviour 
of the Bloch vector is thrown into disorder. In this paper, taking a different viewpoint 
compared with ones that we have been used to, we investigate quasiperiodicity of the 
Bloch vector's trajectories. Using the concept of the quasiperiodicity, we can explain the 
confused behaviour of the system as an intermediate state between periodic and chaotic 
motions. The revealed quasiperiodic structures hidden in the thermal fluctuation of the 
Bloch vector have close relations with topics of the number theory, which is one of the 
oldest branches of the pure mathematics. This appearance of unexpected natures of the 
thermal JCM itself is our motivation for this work. 

In the current paper, we discuss the following three facts: 

• We draw an analogy between the dynamics of the Bloch vector and that of a com- 
pressible fluid. We show that the time evolution of the Bloch vector is equivalent 
to a compressible inviscid flow that has no vorticity. 

• If we adjust the time displacement At properly, figures consisting of plotted dots 
at the constant interval At acquire scale invariance under replacement of At by 
sAt, where s(> 1) is an arbitrary real but not transcendental number. This fact is 
derived by making use of the Weyl criterion, which gives a necessary and sufficient 
condition for a sequence of real numbers to be uniformly distributed modulo unity. 

• We can compute values of the time variable t, which let (^(t)) (the absolute value 
of the z-component of the Bloch vector) be very small, with the Diophantine ap- 
proximation (a rational approximation of an irrational number). 

In Ref. [18], Azuma examines a histogram of {S z (nAt) : n — 0, 1, N}, where S(t) is 
the Bloch vector of the thermal JCM. Plotting the variance of the histogram of samples 
{S z (nAt)} against the inverse of the temperature /3[= 1/(/cbT)], Azuma finds a scaling 
property. In the current paper, we also take a constant interval At for the time variable 
t, and this prescription takes an important role. Thus, the present paper is a sequel of 
Ref. [18]. 

This paper is organized as follows: In Sec. 2, we give a brief review of the thermal 
JCM. In Sec. 3, we examine trajectories of the Bloch vector numerically. In Sec. 4, we 
explain the quasiperiodicity, which we can observe in the trajectories of the Bloch vector. 
In Sec. 5, we draw an analogy between the Bloch vector and the compressible fluid. In 
Sec. 6, we discuss the scale invariance of the figures generated as discrete plots of the 
Bloch vector. In Sees. 7 and 8, we investigate a graph of S x (t) versus the inverse of the 
temperature /3 for the time t such that (^(t)] <C 1 by the numerical experiments and the 
perturbative evaluation, respectively. [S(t) stands for the Bloch vector.] In Sec. 9, we give 
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brief discussions. In Appendix A, we give physical meanings of quantities of a fictitious 
fluid, for instance, the density p, the pressure p and the external force per unit mass K, 
which we introduce for the analogy between the dynamics of the Bloch vector and that of 
the compressible fluid. In Appendix B, we consider how to build the Hamiltonian, which 
yields dynamics of the compressible fluid with zero vorticity. Moreover, we consider the 
reason why the trajectory of the Bloch vector appears to intersect itself as shown in Sec. 3. 
In Appendix C, we give a proof of the Weyl criterion. In Appendix D, we consider physical 
transient spectra of the atom in the cavity. In Appendix E, we prove two theorems, which 
are related to the Diophantine approximation. 



2 A brief review of the thermal JCM 

In this section, we give a brief review of the thermal JCM. To describe its dynamics, we 
use the notation of Ref. [18]. 

The Hamiltonian of the JCM is expressed in the form, 

H = — ujqCi z + hua)a + hg(a + a + a-cJ), (1) 
a± = -(a x ±ia y ), (2) 

^=(10)' o 1 )' ff * = (o -1)' (3) 

[a,a f ] = l, [a, a] = [a f , a f ] = 0, (4) 

where cr x , a y and a z given by Eqs. (2) and (3) are the Pauli matrices acting on atomic 
state vectors, and at and a given by Eq. (4) are the photon creation and annihilation 
operators, respectively. We write the state of the two-level atom as a two-component 
column vector. We describe the state of the cavity field as a superposition of number 
states of photons. 

Let us divide the JCM Hamiltonian defined in Eq. (1) into two parts as follows: 

H = k(C 1 + C 2 ), 
Ci = u(-a z + a 1 a), 

C 2 = g(a+a + a_a j ) - ^-a z , (5) 

where Auj = u — ujq. Then, we obtain the following relation: 

[Ci,C 2 ] = 0. (6) 

Because we can diagonalize the Hermitian operator C\ at ease, we take the following 
interaction picture for describing the state of both the atom and the field. First, we write 
a state vector of the whole system in the Schrodinger and interaction pictures as IV's(^)) 
and \ijji(t)), respectively. Second, assuming \ipi(0)) = |-?/>s(0)), we define \ipi(t)) as follows: 

|^(*)> = exp(iCi*)|^s(*)>- (7) 
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Thus, because of Eq. (6), the time evolution of \ipi{t)) is given by 

Mt)) = umm, (8) 

where 

U(t)=exp(-iC 2 t). (9) 
Here, we give an explicit form of U (t) as follows: 

m = e XPh ,(-y A r /2 )] 

~ (-l)"t 2n / (D + g 2 ) n \ 
= h (2n)! I D»J 

+ g (_i)(_l)n t 2n+l ^ _(A^/2)( J D + ^ 2 ) n (?aL> n ^ 



^ (2n + l)! V ga\D + g 2 ) n (Au/2)D 



uqo = cos(tJ D + g 2 ) + -Aw 



uw = -iga 



2 ' 
sin(tv^D) 



(10) 



D = (^ff + g^a, (11) 



i. sin(tvADT^) 



t sm(V£Tg) 

VD~+~g z ' 

/ /ttn « * shift v/ZT) . . 

cos tV£> - -Aw — V J_ ; . 12 
2 v^D 



We put the whole system into the following initial states: 

Pap(0) = p A (0)®p P , (13) 

Pa(0)= £ P A,iM\i)AA{j\, (14) 
i,je{o,i} 

exp(— Bftuia} a) 

pp = 



Trexp(—(3hua^a) 
= (1 -e- /3 ^)exp(-/3Wa), (15) 

|0>a = ( o ) ' I 1 )a=(?)> (16) 

where {pA,ij} represent an arbitrary 2x2 Hermitian matrix with nonnegative eigenvalues 
and trace unity. The indices A and P stand for the atom and the photon, respectively. The 
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density operator pp given by Eq. (15) represents the thermal equilibrium state weighted 
by the Bose-Einstein distribution for the inverse of the temperature /3[= l/(k^T)]. 

After these preparations, the time evolution of the atomic density operator is given by 



PA(t) = PA,ij(t)\i) AA(j\, 

i,je{o,i} 



(17) 



PA,ij(t) = PA,ki(°) A ki,ij(t) for i,j E {0, 1}, 

fc,Ze{0,l} 



(18) 



A k i,ij{t) = A (i|Tr P [[/(fj(|£;) AA (Z| <g> pp)U*(t)]\j) A . (19) 

Here, we pay attention to the following facts. Clearly, relations PA,io(t) — Pa,oi(^)* an d 
Pa,ii(^) = 1 — PA,oo(t) hold. Thus, we need to estimate only two components, Pa.ooW an d 
Pa,oi(0> so that 



^4oo,oo(^) 
^n,oo(*) 

4)1,01 (*) 



(1 _ e -,r n g (Ac^/2) 2 + g\n + 1) cos 2 (t y/D(n + 1)) ^ 



n=0 



n 



(1 _ e -^) Y\cos{tyl D{n + 1)) + ^-Alj 



n=l 

oo 



sin 2 (t^H) ^_ n0%w 
D(n) 



D(n + 1) 

/3hu 



sin(ty/D(n+l)). 



n=0 



D(n + 1) 



r~-~ % y 

x[cos(ty D(n)) + -Au- 



sin(tyjD(n)). 



4>i,oo(f) = Ai O ,oo(0 = 4)o,oi (t) = Aio,oi (t) = An,oi(*) = 0, 
D(n) = (^) 2 + ,V 



(20) 

(21) 
(22) 



Then, we introduce the Bloch vector S(t) = (S x (t), S y (t), S z (t)), which provides us a 
visual description of the dynamics of the atomic state in a convenient way, 



PA(t)= 1 -[I+S(t).CT], 



S(t) 



( LU(t) Lt(t) \ 
V L<&(t) J 



5(0) + 



/ o \ 



V L<£L(t) J 



(23) 
(24) 



L$,(t) = Re[A 01 , i(t)], 

L ( H(t) = Im[A i,oi(t)], 

^AuW = 4)0,00 (t) — 4l,00 (t), 

L ( Z(t) = 4»,oo(t) + 4i,oo(0 - 1. 



(25) 
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3 Trajectories of the Bloch vector 



In this section, we examine the time evolution of the Bloch vector given by Eqs. (24) and 
(25) numerically. From now on, for simplicity, we assume Au = 0. Then, we can rewrite 
the equations that represent the time evolution of the Bloch vector as follows: 



(26) 





( L^t) 






( ° \ 


S(t) = 


o MO 





5(0) + 







V o o 


Ls(t) J 




V Ht) J 



Ls(t) 
U{t) 



(1 - e~ fj ) J2 cos(V^TTt) cos( v / ^t)e- n/3 , 

n=0 

1 e 2 " - 1 ^ 



^ cos(2v^t)e- n/3 , 



n=l 

\2 oo 



-1(1 - e-P) + £ cos(2^)e-™^, 

Ze n=l 



(27) 



where we assume w 7^ and g ^ 0, and we replace parameters /37iu; and with /3 and 
t, respectively. They imply that the time t is in units of | ^7 1 1 and the inverse of the 
temperature /3 is in units of (hu)^ 1 . As a result of these replacements, the Bloch vector 
S(t) depends only on two dimensionless variables, t and (3. 

Here, putting the initial state of the atom into (1/v / 2)(|0)a + |1)a), we obtain the 
initial Bloch vector 5(0) = (1,0,0) and its time evolution 



S(t) 




(28) 



Equation (28) tells us that the Bloch vector S(t) always lies on the xz-plane Vt(> 0), so 
that it is convenient for tracing the trajectory of S(t) as time passes. Thus, from now on, 
we only examine the case where the Bloch vector is given by Eqs. (27) and (28). 

Figure 1 shows a trajectory of the Bloch vector S(t) given by Eqs. (27) and (28) for 
j3 = 1.0 and < t < 250. To obtain Fig. 1, we replace an infinite summation J2^=o with 
a finite summation X^=o f° r an actual numerical calculation of Li(t) defined in Eq. (27). 
For calculating L±(t) numerically, we give a similar treatment. Throughout this paper, 
whenever we carry out numerical calculations of Li(t) and L^t), we always apply this 
approximation to them. At the end of this section, we argue numerical errors given rise 
to by this approximation. 

Looking at Fig. 1, we feel that both a norm and a direction of S(t) change hard at 
random and it is in a state of disorder. With a careful observation of the trajectory, we 
notice that the figure drawn by S(t) is nearly symmetrical with respect to the z-axis and 
the trajectory lies within a particular limited area on the xz-plane. However, it is difficult 
for us to find a regular form and order any more from Fig. 1. 

As shown in Fig. 1, introducing the non-zero temperature and its thermal fluctuation, 
we can observe that the behaviour of S(t) becomes in disorder. However, if we look at 
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Figure 1: A trajectory of S(t) given by Eqs. (27) and (28) for < t < 250 and (3 = 1.0. 

the time evolution of S(t) from a new viewpoint, which is different from ones that we 
have been used to, we can find novel structures hidden in the trajectory of S(t). We show 
these secret structures in the following paragraphs. 

We plot trajectories of the Bloch vector given by Eqs. (27) and (28) on the res-plane 
as S(t) = (S x (t), S z (t)) in Figs. 2, 3 and 4. In these figures, we plot S(t) at constant time 
intervals, so that the time variable takes discrete values as t n = nAt for n = 0,l,2,...,iV. 
Moreover, we have the number of dots plotted iV as large as possible. Turning our eyes 
on Figs. 2, 3 and 4, we understand that the trajectories of the Bloch vector reveal regular 
forms, distinct orders and novel structures, that we have not been able to discover before. 

Magnifying the upper right part of Fig. 3, we obtain Fig. 5. And magnifying the 
upper right part of Fig. 5, we obtain Fig. 6. Looking at Figs. 2, 3, 4, 5 and 6 carefully, 
we notice that points of \S Z \ <C 1 increases in number as /3 becomes smaller, that is, the 
temperature becomes larger. Moreover, we observe the following facts: On the one hand, 
points of \S Z \ <C 1 for f3 3> 1 satisfy \S X \ ~ 1. On the other hand, points of j<S^| 1 for 
< f3 <C 1 are distributed widely in the range of — 1 < S x < 1. In Fig. 7, we examine 
this tendency of the Bloch vector S(t) numerically. 

How to draw a graph of Fig. 7 is as follows: First, we consider a set of samples, 

{S 2 (t n ) :ne {0,1, 2, ...,N}}, (29) 

where to,ti, ...,t]y are discrete values of the time variable generated with the constant 
interval At as mentioned before. Second, we take a certain positive number being small 
enough, < e « 1. Third, we select values of t n , each of which satisfies \S z (t n )\ < e. 
Because of Eqs. (27) and (28), we can regard S x (= L\) as a function depending only on 
two variables, t n and (3. Fourth, we plot (f3,S x (t n )) such that ^(i^)! < e and obtain 
the graph of Fig. 7. (In Sec. 7, we explain details for drawing Fig. 7, for example, how 
to choose e and iV properly for given ft such that 0.5 < /3 < 5.0, with giving a concrete 
example.) Figure 7 shows dependence of S x on j3 on condition that \S Z \ 1 holds. In the 
graph of Fig. 7, we can observe successive splits leading to different branches. It gives us 
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Figure 2: A figure consisting of dots of S(t) given by Eqs. (27) and (28). We plot them 
at the constant time interval At = 3.5 for j3 = 2.0. The number of dots are equal to 
N = 128 000. We give all points a diameter being 1/1000 of the width of the whole 
graph. 



S z 




-0.7 



Figure 3: A figure consisting of dots of S(t) given by Eqs. (27) and (28). We plot them 
at the constant time interval At = 3.5 for (3 = 1.0. The number of dots are equal to 
iV = 384 000. We give all points a diameter being 1/1000 of the width of the whole 
graph. 
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Figure 4: A figure consisting of dots of S(t) given by Eqs. (27) and (28). We plot them 
at the constant time interval At = 3.5 for /3 = 0.5. The number of dots are equal to 
N — 1 024 000. We give all points a diameter being 1/1000 of the width of the whole 
graph. 




Figure 5: A figure consisting of dots of S(t) given by Eqs. (27) and (28). We plot them 
at the constant time interval At = 3.5 for (3 — 1.0 and iV = 12 800 000. The horizontal 
and vertical ranges are given by 0.7 < S x < 1.0 and —0.1 < S z < 0.0, respectively. This 
figure corresponds to an enlargement of the upper right part of Fig. 3. The number of 
dots plotted actually for this figure is equal to 274 692. We give all points a diameter 
being 1/1000 of the width of the whole graph. 
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Figure 6: A figure consisting of dots of S(t) given by Eqs. (27) and (28). We plot them at 
the constant time interval At = 3.5 for j3 — 1.0 and N — 1 024 000 000. The horizontal 
and vertical ranges are given by 0.95 < S x < 1.0 and —0.02 < S z < 0.0, respectively. 
This figure corresponds to an enlargement of the upper right part of Fig. 5. The number 
of dots plotted actually for this figure is equal to 507 596. We give all points a diameter 
being 1/1000 of the width of the whole graph. 




Figure 7: Taking a certain positive number being small enough as < e C 1, we plot 
(f3, S x (t n )) for discrete time variable t n = nAt for At = 3.5 on condition that ^(in)] < e 
holds. The horizontal axis represents (5 for the range of 0.5 < /3 < 5.0 with the logarithmic 
scale. We plot all points at a constant interval 0.025 for /3. We give them a diameter being 
4/1000 of the width of the whole graph. (We explain how to draw this graph precisely in 
Sec. 7.) 
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an impression that the Bloch vector S (t) hides further secret orders and structures inside 
itself. 

Closing this section, let us argue the numerical precisions of Li(t) and L 4 (t). As 
mentioned before, whenever we carry out the numerical calculations of Li(t) and L 4 (t) 
defined in Eq. (27), we replace the infinite summation Y^=o with the finite summation 
Z)n=o> which contains the first 151 terms. Here, we estimate numerical errors caused by 
this treatment. For example, we can evaluate the upper bound of numerical errors for 
Li(t) as 



< 


(1-e"' 


oo 

3 )l E < 






n=151 






oo 


< 


(1 - e"' 


3 ) E e 






n=151 




(1 - e"' 


3 )e -151/3 




e -151/3_ 





-n/3 



E^ n/3 

n=0 

(30) 



Thus, if we assume 0.5 < (3 < 5.0, we obtain 

£[Li] < 1.63 x 10~ 33 . (31) 

Therefore, our calculation of serves a precision of 33 significant decimal digits to us. 
Similar things hold for L^(t). 

To obtain numerical data for drawing graphs of Figs. 1, 2, 3, 4, 5, 6 and 7, we utilize 
the Fortran compiler, which can manipulate the data type of the real quadruple precision. 
This implies that the compiler sorts out numerical calculations with a precision of about 
33 significant decimal digits, so that it is fit to suppress numerical errors as Eqs. (30) and 
(31). 

To obtain numerical data for drawing graphs of Figs. 13, 14, 15, 16, 17, 22, 23, 
24, 25, 26, 27 and 28, and to carry out calculations required in Sec. 8, we utilize a 
high-performance computer algebra system, the Mathematica. The Mathematica has an 
option allowing us to specify any precision we want, so that we let it execute numerical 
calculations with 50 significant decimal digits for keeping extra enough amount of digits. 
In Sec. 8, we often have to manipulate integers, whose number of decimal digits is equal 
to 15, in computations. Under these circumstances, we need 15 digits for integers and 33 
digits for the decimal representation of the floating point, so that we have to prepare 48 
digits in total for computations over real numbers. Thus, utilizing the Mathematica with 
the 50 significant decimal digits, we can keep away from danger caused by the roundoff 
and truncation errors. Therefore, we can obtain numerical values of Li(t) and L±(t) with 
33 significant decimal digits at least, and we can obtain stable results of the numerical 
calculations. 
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4 Quasiperiodicity in the dynamics of the Bloch vec- 
tor 



In this section, we discuss the quasiperiodicity of the Bloch vector, which is given by 
Eqs. (26) and (27). Here, we examine the behaviour of the Bloch vector under the low 
temperature limit. Assuming that the temperature is low enough as (3 = l/ik^T) ^> 1 and 
neglecting the terms of order 0((e~ /3 ) 2 ), we obtain approximations of L^t) for i — 1,3,4 
as follows: 

Li(f) ~ (1 -e"' 3 ) cost + e~ /3 cos(v / 2Y) cost, (32) 
L 3 (t) ~ I(l- e -/») + Icos(2t), (33) 
Ut) ~ -1(1 - e -/») + 1(1- 26-^)008(2*). (34) 

The above approximations of L 3 (t) and L±{t) are periodic functions of t with period n. 
In contrast, the function Li(t) is made of cost and cos(v / 2t), that is, a complex structure of 
trigonometric functions with periods 2% and V2n, so that its behaviour is unpredictable. 
Thus, from now on, we concentrate on analyzing the function L\{t). Equation (32), the 
approximation of Li(t), has two angular frequencies, oui — 1 and w 2 = v^2, so that the 
ratio of w 2 to cui is irrational as u; 2 /Vl = V2- We call m angular frequencies wi, o> 2 , to m 
are incommensurate if no one of the angular frequencies Ui can be expressed as a linear 
combination of the others using coefficients that are rational numbers. And if the system 
has incommensurate angular frequencies, we call it quasiper iodic. 

The notion of the quasiperiodicity means an intermediate state between periodic and 
chaotic ones [19, 20, 21, 22, 23, 24]. To understand the quasiperiodicity of Eq. (32) giving 
the approximation of Li(t), we examine a dynamical system that has two incommensurate 
angular frequencies, ui\ and u 2 . In the case of two- angular-frequency quasiperiodic motion, 
we can describe a dynamical variable /(t) as a function of two independent variables, 
G(ti,t2), for example, 

f(t) = Li (t) = (1 - e-1 3 ) cos t + e _/3 cos( y/2t) cos t, (35) 

G(ti,t 2 ) = (1 -e^costi + e~ /3 cos(v / 2t 2 )costi, (36) 



G(h + 2vr, t 2 ) = G(h, t 2 + V2n) = G(h, t 2 ), (37) 
f(t) = G(t 1 ,t 2 )\ ti=t2=t . (38) 
Letting u\ — 1 and w 2 = y/2, the equation of the form, 

m\UJ\ + m 2 w 2 = 0, (39) 

does not hold for arbitrary integers, m\ and m 2 , except for m 1 = m 2 = 0, so that u)\ and 
u>2 are incommensurate. Thus, the dynamics of the system is specified by two independent 
angles 61 and 9 2 as 

G{W (40) 
OJi 0J2 
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where < 9i < 2n for i — 1,2. This implies that the orbit of the motion lies on the torus 
embedded in the phase space. From Eq. (38), we understand that f(t) is the quantity 
obtained along the orbit of Q\ = u\t and 6 2 = uj 2 t being realized on the torus in the phase 
space. In the following, we consider these matters in a general and precise manner. 

Let us consider an arbitrary Hamiltonian H(qi,q 2 ,pi,p2,t), whose number of degrees 
of freedom is equal to two. Because the system evolves according to this Hamiltonian, 
the trajectory that the system follows is described with Hamilton's equations, 

dH({q,},{p,},t) 

OPi 

= jjHMAiAA fori = i, 2 . (4i) 

dqi 

The phase space of the system is four-dimensional and its co-ordinate system is given 
by (<?i, <?2,Pi,P2)- Then, we remember the Poincare-Cartan integral invariant and Liou- 
ville's theorem, which are fundamental results in the field of the analytical mechanics 
[21, 24, 25]. They guarantee that the density of the system points in an infinitesimal 
volume element is preserved while they are traveling with time according to the canonical 
equations of motion. In other words, the volume of a region of the phase space is invari- 
ant under evolution with Hamilton's equations of the motion in time. Thus, for example, 
dqi A dpi, dq 2 A dp 2 and (dqi A dpi) A (dq 2 A dp 2 ) are preserved along the trajectories in 
the four-dimensional phase space. 

Now, we assume that the motion of the system is bounded in the four- dimensional 
phase space. Furthermore, we suppose the Hamiltonian system to be integrable. Thus, 
the system has two functions Ii({qj}, {pj}) for i = 1,2, which are in involution and 
independent. Moreover, these functions are time-independent and preserved. Hence, 
they satisfy the following commutation relations: 

[J 1 ,/ 2 ] = 0, [h,H] = for i = 1,2. (42) 

In the above commutation relations, the bilinear operation [u, v] denotes the Poisson 
bracket, so that ij for i — 1, 2 are constants of the motion. 

Here, let us define ij for % — 1, 2 at initial time (t = 0) as follows: 

^(0) = li({ qj (0)}, {Pj(0)}) = constant for % = 1, 2. (43) 

This implies that the motion of the system is confined to a two-dimensional surface em- 
bedded in the four-dimensional phase space (qi,q2,Pi,P2)- In addition, the surface is 
determined by the two integral of the motion, ii(0) and ^(0), and all possible states of 
the system have to lie on it. 

In this section, we often use two technical words, a constant of motion and an integral 
of motion. The constant of motion is a quantity that is conserved through the motion. On 
the other hand, the integral of motion is a function which is constant along a trajectory in 
the phase space. Clearly, an integral of the motion is a constant of the motion. However, 
a constant of the motion is not necessarily an integral of the motion. 

From now on, we let -Mj(o) represent the above surface specified with ii(0) and ^(O). 
Then, we suppose /iq G -Mj(o)to be an initial point of the system at t = in the phase 



14 




Figure 8: The sequence of the open neighbourhoods, cr^e), a^ A t), C//(2At), cr^NAt), 
which lies on the surface M.i(q). Because an area of Aii(o) is finite, there has to exist an 
integer JV(> 0) such as (7 w (e) n cr^NAt) ^ <p. 

space. Moreover, we define an open neighbourhood with centre /x an d radius < e < 1 
on the surface A4/( ), such that 

°Mo( e ) = {A* g Mi(q) : d(n,fi Q ) < e}, (44) 

where d(n, //q) represents a distance between \x and //q on M.i(q). 

A point in the phase space takes a set of values for the canonical co-ordinates and 
momenta. It evolves through the phase space in time. We let the point travel from //o 
at t — to //(At) for a period of time At, where < At < 1. Furthermore, we let 
°>(At) be a set of points which come from <r w (e), the open neighbourhood of // , along 
the trajectories of the system evolving for a period of time At. Clearly, //(At) G c M (At) 
and cr^ At) e .M/(o) hold. 

Because of the Poincare-Cartan integral invariant and Liouville's theorem, an area 
of cr^At) has to be equal to that of o" M() (e) on the surface Aii(p). Then, we repeat this 
time evolution many times with a constant interval At. At the same time, we adjust 
some parameters as follows. For an arbitrary but fixed time interval At, we can always 
take a sufficiently small radius e(> 0) of the open neighbourhood <7 Mo (e), such as any two 
successive open sets never intersect in the sequence of c Mo (e), c^At), C/i(2At), .... Namely, 
making e(> 0) small enough, we can let the following relations hold, 

o" M (nAt) H <7 M((n+1) At) =<P for n = 0, 1, 2, .... (45) 

From the above considerations, the sequence of the open neighbourhoods, u^nAt) f° r 
ra = 0,l, 2,..., lies on the surface Aii(o) as shown in Fig. 8. 

Here, we remember the following fact. We assume that the motion of the system is 
bounded in the phase space. Thus, an area, where all points of actually possible states of 
the system belong to, is finite. This implies that the sequence of the open neighbourhoods, 
<7 W) (e), <r M (At), cr M (2At), lies on the region of finite size. Hence, there has to exist a finite 
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Figure 9: A two-dimensional co-ordinate system {6i,6 2 ) on the surface Aii(o). 
positive integer N{> 0) such that 

°mo ( e ) n At) 7^ 0- (46) 

This is because every open neighbourhood of cr^nAt) for n = 0, 1, 2, ... has a finite equal 
area, which is determined by radius e. If the sequence of open neighbourhoods covers the 
entire finite region of actually possible states on the surface A4i(o), there must exist a 
certain finite positive integer N that satisfies Eq. (46). 

The results obtained above hold even if we take arbitrary small At and e. Hence, we 
understand that the motion of the integrable system bounded on the surface Mi(o) has 
to be periodic. This conclusion implies that the motion of the system is periodic along 
the trajectory of the time evolution being determined by Hamilton's equations. 

Moreover, we can extend this statement as follows. At first, we assume that ii(0) 
and ^(0) are the integrals of the motion of the system. Next, we regard these inte- 
grals as general angular momenta h{{lj},{Pj}) for i — 1,2. Let Q\ and 6 2 be canonical 
angular variables of these angular momenta, so that {6i,I\) and {6 2 ,I 2 ) form canonical 
co-ordinates. Then, putting parameters along the loop as Q{ G [0, 27r) for % — 1,2, we 
can consider both variables of 6 X and 6 2 to be periodic, 

e t = idit (mod 2tt) for % = 1, 2. (47) 

We can regard Eq. (47) as a definition of o>j for i — 1, 2. In the following paragraphs, we 
explain the reason why Eq. (47) holds more precisely. 

On the surface jMj(o) being determined by the two integrals of the motion, ii(0) and 
J 2 (0), we can put a two-dimensional co-ordinate system {61,62), where 6\ and 6 2 are 
independent of each other, as shown in Fig 9. At this moment, we suppose 61 and 6 2 do 
not need to be periodic. Here, we pay attention to the following fact. Using the Poisson 
bracket, we can describe a differential operator with respect to a variable &i on the surface 
Mi( ) as 

dl d dl d d 

^^-^wr^w (48) 
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This implies that a derivative of a function f(6i, 9 2 , -Tl(O), ^2(0)) in 6\ direction at an 
arbitrary point on the surface Aiz(o) is given by 

D I J = -[I l ,f] for i = 1,2. (49) 

Because we suppose [ii, J 2 ] = in Eq. (42), the above two differential operators commute 
obviously as 

D h D h -D h D h = -[/i,[/ 2 ,-]] + [h,[h,-}} 
= - D [h,h] 

= 0. (50) 

Moreover, D Tl and D l2 are linearly independent because Q\ and 6 2 are independent 
variables. We explain this fact in detail as follows. Let us suppose there exist non-zero 
coefficients a± and a 2 such as aiDj 1 + a 2 Di 2 = 0. Furthermore, we consider an arbitrary 
function f(9i, 6 2 ) defined on .Mj(o), and we obtain an identity, (aiD^ + a 2 Di 2 )f(9i, 6 2 ) = 
0. This implies 

Mi + a 2 I 2 , f] = (<*iQf + a2 W )f(0i, 9 2 ) = (51) 

The above identity has to hold for any function f{6i,6 2 ) defined on A4i(p), so that we 
obtain 0:1 = 0:2 = 0. Hence, we can conclude that Di 1 and Dj 2 are linearly independent. 

Namely, {D Tl ,D l2 } is a set of linearly independent and commuting differential oper- 
ators. Moreover, we can regard these differential operators as generators of the Abelian 
group of displacements on Aii(p). Here, we put a reference point (an origin) of the co- 
ordinate system (6*1,6*2) at an initial point /xo for t — on Aii(p). Moreover, we consider 
(t±Ii + t 2 I 2 ) to be a generator of a canonical transformation. This canonical transfor- 
mation causes the following displacement of the reference point. We let fx Q denote the 
two-dimensional generalized co-ordinates of fj, G A^j(o), and fi(0) = fi(9i,9 2 ) denote 
the two-dimensional generalized co-ordinates specified with (6*1, 6 2 ) on A^j(o)- Then, we 
obtain 

Mo = A*l«=(o,o) — ► ^(Mo) = Mo- E M 7 i»A*(0)] 

ie{i,2} 

/i(o,o)+ E *■ 

ie{l,2} 



80^ 



(52) 

e=(o,o) 



From the above equation, we notice that derivatives concerning the general angular mo- 
menta U for i — 1,2 are equal to zero during the displacement of g t (fJ>o) from /i . Thus, 
<?*(/io) is always confined on the surface M.i(o). 

Here, we pay attention to the following fact. A change of an arbitrary function F with 
time t is given by the following Hamilton's equation: 

dF dF 

If the function F is independent of time t in an explicit manner, we obtain 

f = [F lB ]. (54) 
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Figure 10: We can consider the surface M.i(o) to be topologically equivalent to a two- 
dimensional torus, on which every point is specified with the co-ordinate system (8i,8 2 )- 

The above equation implies 



We regard the above equation as a canonical transformation whose generator is given by 



another new Hamiltonian. Then, regarding Eq. (56) as its time evolution, we can apply the 
Poincare-Cartan integral invariant and Liouville's theorem to the canonical transformation 
caused by £ie{i,2} Uk- 

Thus, when the open neighbourhood with centre /x and radius e is mapped to an 
open neighbourhood with centre g t (f J Lo), its area is preserved along a path of the canonical 
transformation. This means that there exists an certain positive integer N, such as an 
open neighbourhood with centre (g t ) (2 ' N (^o) intersects the initial open neighbourhood 
with centre /x in the sequence of fx — > <?*(/io) — > <7*(<7*(a*o)) — > ■■■ —> (g*)® 1 * (l^o) ■ 

In this discussion, we can suppose t = (£1,^2) to be an arbitrary two-dimensional real 
vector. Hence, we understand that the two angles 81 and 82 are independent and periodic 
in times t\ and ti, respectively. Thus, we obtain $i = for % — 1,2. Then, reminding 
Eqs. (35), (36), (37) and (38), we set t\ — t 2 — t and obtain Eq. (47). From the above 
discussion, we can consider that the surface .Mj(o) embedded in the four-dimensional 
phase space to be topologically equivalent to a two-dimensional torus, T 2 = S 1 x S 1 , as 
shown in Fig. 10. 

Here, let us think about the special case 002/001 = 1/2. Then, an orbit winds around 
the torus twice in the short direction 8\ and simultaneously once in the long direction 8 2 
as time progresses. Thus, the orbit closes on itself after finite period of time as shown in 
Fig. 11. Hence, actually possible states of the system are confined on the one-dimensional 
closed path lying on the torus. 

However, if we think about the case where uj^I^x is irrational, for example, oj 2 /uji = 
y/2, the situation with respect to the orbit of the system becomes different drastically. 
If U2/0J1 is irrational, the orbit never closes on itself as shown in Fig. 12. As time goes 




(55) 



dtH. Then, a canonical transformation of /j,q caused by a generator J2ie{i,2} can be 
written as 
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Figure 11: In the case of uj 2 /uj 1 = 1/2, an orbit moves around the torus twice in the short 
direction 9i and once in the long direction 9 2 as time proceeds, so that the orbit forms a 
finite closed path. 




Figure 12: In the case of u 2 /ui = V2, an orbit never closes on itself in a finite time 
period. The trajectory becomes dense on the torus as time proceeds to infinity. 
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to infinity, the orbit on the torus will eventually come arbitrary close to every point on 
the toroidal surface. Thus, the orbit fills up the torus in the phase space finally, and the 
trajectory of the system becomes dense and uniformly distributed on the torus. 

Furthermore, if U2/001 is irrational, the trajectory of the system has both fractal and 
ergodic properties that periodic systems never show. We explain these characteristics of 
quasiperiodic systems in the following. At first, let us consider an arbitrary continuous, 
smooth and bounded function defined on the surface Mi(o), such as /(//) V/i G M.i(o). 
Next, we define the following two averages of /(//). 

• The time average of /(//). Getting a reference point as a starting-point (an initial 
point of t — 0), we take an average of /(//) along the path, which is parameterized 
with time t G [0, +00) as 

(f)M = lim I f T f(wt)dt, (57) 

T^+oo 1 Jo 

where to = (001,002). In the above calculation, we put the reference point /i at the 
origin of the two-dimensional co-ordinate system 6 = (9i,9 2 ) = (0,0). 

• The average of /(//) over the phase space. We take the average of /(//) over all 
points of /j, belonging to the surface -Mj( ) as 

(/> "'« = (2^i dh L dhm - (58) 

Then, we can show the following fact. If U2/001 is irrational, the above two averages 
become equivalent to each other [24], 

(f)t = (f)M m - (59) 

The reason why is as follows. 

We introduce a Fourier-series expansion of f(0) on M.i(o) as 

f(o)= E hz lke i (6°) 

k=(k u k2)ez 2 

where Z represents the set of all integers. Taking care of 

— / d6e m = 4,o, (61) 

27T JO 

we obtain 

(f) Mm = /o, (62) 
at ease. On the other hand, we can rewrite the time average as 



(f)M = ^^[dt E E 

1 " ' fcie{o,±i,±2,...}fc 2 e{o,±i,±2,...} 



= fo + y fk hm - f dt < ,A - ; 

k7o J t ^tJ 



1 e ik-uiT _ 1 
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In the derivation of the above equation, we use the fact that k ■ uj ^ Vfc 7^ if U2/0J1 is 
irrational. Moreover, using the following formula, 

lim - — = forfc-w^O, (64) 

we arrive at 

</>t(^o) = fo- (65) 

From Eqs. (62) and (65), we can conclude Eq. (59). 

These results tell us the following. If u 2 /oji is irrational, the trajectory of the system on 
Aii(o) has the Hausdorff dimension Dh = 2, which is related with the integral of a function 
defined on the surface Adito)- However, the topological dimension of the trajectory on 
A4i(o) is obviously given by = 1. Thus, we obtain < Dh and we have to conclude 
that the trajectory on -M/(o) is fractal [26, 27]. Moreover, the relation (f) t (fi ) = (f)M I{0) 
implies that the time average and the average over phase space are equivalent to each 
other, so that the trajectory has ergodicity. From these discussions, we can expect that 
the trajectory drawn by the Bloch vector has fractal properties if u^/u^i is irrational. We 
examine this point in Sec. 6. 

Here, we pay attention to the following fact. As mentioned before, if U2/0J1 is irrational, 
the trajectory never closes on itself in the phase space. However, taking the long-term 
limit, the trajectory on the torus in the phase space will eventually comes arbitrary close 
to every point on the toroidal surface. Hence, if we take an arbitrary small positive 
number e(> 0), there actually exists time t such as 

d(fjLo, n(t)) <e. (66) 

The question is how to estimate time t, which satisfies Eq. (66). Let us think about 
obtaining a rational approximate number of irrational w 2 /wi by its continued fraction 
expression. For example, in case of w 2 /^i = V%, we can describe its continued fraction as 

V2 = l + ^ . (67) 

2 + - 

2 + — 

By truncating the above expression with the first some terms, we can obtain the rational 
approximate number p/q with high precision such as 



(This fact is assured by the theorem of Roth 

Thus, because of U2/0J1 = \/2 ~ p/q, we obtain pco 1 ~ qco 2 - This implies the following. 
If the orbit winds around the torus pth times in the short direction 0\ and simultaneously 
gth times in the long direction 62, it arrives at a point which is very close to the initial 
point /io- We investigate this question in Sec. 8. 

In this section, so far, we suppose the low-temperature limit and evaluate the con- 
tribution to Li(t) up to 0(e^ 13 ). Here, we let the temperature become higher slightly 
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and estimate the contribution to L\(i) up to 0((e~") 2 ). Then, we can regard L\{t) as a 
superposition of cost, cos(v^2£) and cos(y / 3£), and it causes quasiperiodicity of three in- 
commensurate angular frequencies. Moreover, if we make the temperature become much 
higher and estimate the contribution to Li(i) up to 0((e _/3 ) 3 ), Li(t) has four angular 
frequencies, u>i — 1, ui 2 — \/2, UJ3 = V% and U4 = 2. Then, cj 4 /wi is equal to two, so that 
it is rational. Let us think about the four-dimensional torus T A = S 1 x S 1 x S 1 x S 1 , 
where the orbit of the motion lies, and its two-dimensional submanifold (61,64). Then, 
the orbit is periodic on the submanifold (61,64). So that, as the temperature becomes 
higher, the motion of the Bloch vector causes a complicated trajectory, which is a mixture 
of periodic and quasiperiodic orbits. 

Finally, we note the fact that L\(t) has countably infinite incommensurate angular 
frequencies, u n = y/n for n = 1,2,3,5,.... This means that the system has countably 
infinite angular momenta as integrals of the motion, I n for n = 1,2,3,5,.... Thus, we 
have to conclude that we cannot describe the motion of the Bloch vector as the n-body 
problem in classical mechanics, which has a finite number of degrees of freedom. Hence, 
we can expect the motion of the Bloch vector to be equivalent to a classical field theory, 
which has a countably infinite number of degrees of freedom. In fact, we confirm that the 
motion of the Bloch vector is equivalent to the dynamics of a compressible fluid in Sec. 5. 

5 Equivalence between the Bloch vector and com- 
pressible fluid 

In this section, we discuss equivalence between the dynamics of the Bloch vector and a 
compressible fluid flow. In Sec. 4, we emphasize that the quasiperiodic motion of the 
Bloch vector has countably infinite angular momenta, so that we cannot describe it as 
the n-body problem of the classical mechanics. 

We derive the Bloch vector S(t) given by Eqs. (26) and (27) from the Schrodinger 
equation, whose Hamiltonian is defined in Eqs. (1), (2), (3) and (4), originally. As ex- 
plained in Sec. 2, for obtaining S(t), what we apply to the quantum mechanical state in 
fact are only providing the initial state of the cavity field at t = as a density matrix of 
a mixed state and tracing out its Hilbert space V£ > 0. 

Because we trace out the Hilbert space of the cavity field, the entanglement between 
the single atom and the cavity field is removed. Thus, the dynamics of the Bloch vector 
has to be described with classical theory. Hence, we can expect that the trajectory of the 
Bloch vector to be described with a classical field theory, which has a countably infinite 
number of degrees of freedom. In this section, we demonstrate that it is exactly equivalent 
to hydrodynamics of a compressible fluid, which is inviscid and has zero vorticity [29, 30]. 

From now on, we describe S(t) as x(t) = (x(t),y(t), z(t)). Thus, we rewrite Eq. (26) 

as 

/ x(t) \ I Li(t)a;(0) \ 

y(t) = £i(f)?/(o) , (69) 

\z(t) J \ L 3 (t)z(0) + L A (t) J 

where |as(0) j < I. We consider Eq. (69) to be the Lagrangian forms of the dynamical 
equations of a fluid. In Eq. (69), we let £c(0) be the initial co-ordinates of any particle 



x(t) 
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of the fluid and x(t) be its co-ordinates at time t. Hence, Eq. (69) expresses the time 
evolution of the fluid particles being put initially inside a unit circle, that is, |x(0)| < f. 
From Eq. (69), the velocity of the fluid particle at point x(t) at time t is given by 



x(t) 




( Li(i)x(0) 
Li(t)y(0) 
V L 3 (t)z(0) + L 4 (t) 



(70) 



Thus, from Eqs. (27) and (70), calculating an initial velocity at every point inside the 
unit circle explicitly, we obtain 



v(t = 0;x(0),y(0),z(0)) = x(t = 0; x(0), y(0), z(0)) 

= Va;(0) such that |cc(0) j < 1. 



(71) 



To obtain the dynamical equations of the fluid, we rewrite Eq. (69) as the following 
relations, 



L^ty iaity 

Substitution of Eq. (72) into Eq. (70) yields 



x(t) = 



(72) 



[LM/LtfMt) 
[iaW/IaitMt) 
V [L 3 (t)/L 3 (t)][z(t)-L 4 (t)]+L 4 (t) J 



(73) 



Therefore, the velocity of the fluid at point x = (x, y, z) at time t is described as 



v(t,x,y,z) 



[L^/L^x 
[L^/L^jy 
V [L 3 (t)/L 3 (t)][z-L 4 (t)] + L 4 (t) 



(74) 



Here, let us point out some characteristic properties of the fluid obtained by the 
analogy of the Bloch vector. From Eq. (74), we obtain 



V • v 



MM ■ Ls(t) 

%(0 Hty 



(75) 



Thus, we understand that V • v depends only on j3 and t. By numerical calculations, we 
can confirm V • v ^ at almost all times t for < /3 < +oo. Thus, we consider that the 
following relation holds in general, 

V • v j£ 0. (76) 
In addition, from Eq. (74), we can obtain the following equation, at ease: 



V x v = V(t,x). 
Then, let us think about the Eulerian equation of continuity, 

| + v.(H = o. 



(77) 



(78) 
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If we assume that the fluid is incompressible, that is, p = constant, we obtain V • v = 
\/(t,x). However, it contradicts Eq. (76), so that we conclude that the fluid has to be 
compressible. Moreover, Eq. (77) implies that there is no vorticity in any portion of the 
fluid at any time. In general, a concept of turbulence closely relates to randomness and 
complexity of vorticity. Thus, the fluid never shows turbulence motion. 

Let us examine dynamics of the fluid further. From Eqs. (74) and (77), we obtain the 
following velocity-potential $: 

v = V$, (79) 
Then, from Eq. (80), we obtain the following relation: 

Thus, using Eq. (79) and (81), we arrive at 

Av = AV$ = VA$ = 0. (82) 

Equation (82) implies that the viscosity of the fluid has no effects on evolution of an 
actual motion of a fluid element. The reason why is given in the following. 
In general, the equations of motion of a compressible fluid are given by 

^ = % + (v ■ V)« = K - -Vp + ^A,, (83) 
Dt at p p 

which are called the Navier-Stokes equations. In Eq. (83), p, p, K and p represent the 
density of the fluid, the pressure, an external force per unit mass and the viscosity of the 
fluid, respectively. However, because of Eq. (82), we can drop the term including p from 
the Navier-Stokes equations. Thus, we understand that the viscosity p has no effects on 
the motion of the fluid, in fact. 

Putting the above considerations together, we obtain the following observations. If we 
regard the Bloch vector S(t) as a fluid, it is compressible and inviscid. Moreover, it has no 
vorticity. By reason of the fact that the fluid is inviscid and has zero vorticity, its motion 
seems to be simple for us. However, this impression is not true for the compressible fluid. 
Indeed, sometimes we can hardly predict the motion of the compressible fluid, even if it 
does not have the viscosity. 

Figures 13, 14 and 15 show the velocity v(x,y,z,t) given by Eq. (74) on the res-plane 
for \x\ < 1 and z < with (3 — 1.0 at t — 0.5, 1.0 and 1.5, respectively. Thus, we can 
think that Figs. 13, 14 and 15 represent the tangent vectors of S(t) at t = 0.5, 1.0 and 1.5, 
respectively. Looking at these figures, we understand that the velocity of the fluid at any 
point at any time changes hard, so that the particle of the fluid draws a very complicated 
trajectory. 

In Appendix A, we give physical meanings of p, p and K defined in Eqs. (78) and 
(83). In Appendix B, we consider how to build the Hamiltonian, which generates the 
Navier-Stokes equations (83). Moreover, we consider the reason why the trajectory in 
Fig. 1 appears to intersect itself. We explain that intersections never violate existence 
and uniqueness of solutions for the Navier-Stokes equations. 
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Figure 13: This figure represents the velocity v given by Eq. (74) at points of \x\ < 1 
and z < on the xz-pl&ne at time t = 0.5 with (5 = 1.0. The lengths of arrows (vectors) 
represent only just the ratios of their norms compared with each other. 




Figure 14: This figure represents the velocity v given by Eq. (74) at points of \x\ < 1 
and z < on the xz-plane at time t = 1.0 with (3 = 1.0. The lengths of arrows (vectors) 
represent only just the ratios of their norms compared with each other. 
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Figure 15: This figure represents the velocity v given by Eq. (74) at points of \x\ < 1 
and z < on the xz-pl&ne at time t = 1.5 with (3 = 1.0. The lengths of arrows (vectors) 
represent only just the ratios of their norms compared with each other. 

6 The scale invariance of the figures generated by 
the Bloch vector 

In this section, we show that the figures drawn in Figs. 2, 3, 4, 5 and 6 are invariant 
under a scale change in At. We draw these figures by plotting the trajectories of the 
Bloch vector at the constant interval At. This scale invariance is manifestation of the 
fractal structure of the Bloch vector's dynamics, as mentioned in Sec. 4. 

To show the scale invariance of these discrete plots, we use some properties of random 
sequences of real numbers, which are uniformly distributed modulo 2rr. In the former 
half of this section, we give two considerations in preparation for explaining the scale 
invariance of the figures. First, we discuss the physical meaning of plotting the trajectories 
of the Bloch vector at the constant time interval At. We examine differences between 
continuous and discrete plots. Second, we discuss the random sequences that are uniformly 
distributed modulo 2-n. After these preparations, in the latter half of this section, we 
discuss the scale invariance of figures generated by plotting the Bloch vector at the interval 
At. 

First, we consider physical meaning of dividing time into a lattice. With a discrete 
time slice At, we have a chain of time co-ordinates t n = nAt for n = 0, 1, 2, N. Then, 
as shown in Figs. 2, 3, 4, 5 and 6, we plot trajectories of the Bloch vector S(t) as 
{S(t n ) : n E {0, 1, 2, N}}. Through this procedure, we introduce a unit of time At. 
From now on, we regard it as a scale of the time variable. Indeed, assuming At to be a 
real but not transcendental number and At > tt, the graphs of Figs. 2, 3, 4, 5 and 6 are 
invariant under the following change of the scale, 

At ->■ sAt, (84) 

where s is a real but not transcendental number and s > 1. We show this fact in the 
latter half of this section. 

As an alternative plan to Eq. (84), we can consider a change of a scale of the continuous 
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time variable, 

t st (s> 0). (85) 

In this case, we handle the time variable t as a continuous quantity. Throughout the 
current paper, we adopt Eq. (84) rather than Eq. (85). Here, we think about a difference 
between the change of the scale of Eq. (84) and that of Eq. (85). On the one hand, 
if we choose Eq. (84), we have to neglect all events which relate to the quantity of the 
time being less than At. On the other hand, if we select Eq. (85), we have to admit 
an infinitesimal time displacement for the system. Thus, taking Eq. (84), we get rid of 
notions of infinitesimal time evolution. In the latter half of this section, we consider At 
to be in the range between ir and 2ir, that is, n < At < 2n. This treatment generates a 
random sequence t n = nAt (mod 2n) for n = 0, 1, N and lets the graphs of discrete 
dots in Figs. 2, 3, 4, 5 and 6 be invariant under the change of the scale At. 

Moreover, the above treatment produces an unexpected effect that removes a history 
of the sequence {S(t n ) : t n = nAt,n G {0, 1, 2, N}} actually. This is because At is 
rather larger, such as n < At < 2n. If we let the time displacement be small enough, 
the sequence of (S(t ), S(ti), S(t 2 ), ...) approximates to the continuous trajectory of S(t) 
well. However, if we take n < At < 2tt, we cannot find neither trace nor history of 
continuous trajectory of S(t) except for its discrete dots. 

Next, we think about random sequences uniformly distributed in [0,27r). At first, 
we study some properties of sequences of real numbers, which are uniformly distributed 
modulo 2tt (but need not be random). We begin by giving the following set Sn, which 
consists of (N + 1) elements, 

S N = { Xj : < Xj < 2n,j e {0,1, N}}. (86) 

Rewriting the elements of Sn as an ordered list (xq,xi, ...,xn), we can consider it to be 
a real sequence. Let us examine conditions for the set Sn to be uniformly distributed in 
the range of [0, 2ir) under the limit of TV — > +oo. If the set S N is uniformly distributed, 
the following relation holds, 



I N 1 i-2-k 

lim — — -Y J <t>{x j ) = — dx(f)(x), (87) 
jv^+oo 7V + 1 31 2ir Jo y ' y 1 



where <f)(x) is an arbitrary real valued function being bounded, differentiable and contin- 
uous for < x < 2n. 

Here, let us give a Fourier-series expansion of (f>(x) as follows: 



+oo 



m=—oo 



c m = ±- ' dx(t>{x)e- imx . (88) 

Z7T JO 

Then, substitution of Eq. (88) into the left-hand side of Eq. (87) yields 

Y N 1 N +oc 
lim — 6(xi) = lim — c m exniimx A . (89) 
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In a similar way, substitution of Eq. (88) into the right-hand side of Eq. (87) yields 

/ dx ^ x ) / dx £ Cme imx (9Q) 

Thus, putting Eqs. (87), (89) and (90) together, we obtain 

a&o/vZT £ E c ™ ex p(^) = ^-/ n E c ™ eMra ( 91 ) 

' m=— oo j=0 rrt=— oo 



Then, we remember that 0(x) is an arbitrary function in Eq. (87). Hence, we can regard 
{c m : m E {±1,±2, ...}} as arbitrary coefficients in Eq. (91), so that we obtain the 
following relations: 

1 N 1 /-27T 

lim exr>(imXj) = — / dxe imx 

N^+ooN+lf^ HV 11 2tt Jo 



r- 

1 [e imx 2 ^ 



x=0 



27T L «m 

= Vme{±l,±2,...}. (92) 

These relations are called the Weyl criterion [28, 31, 32, 33, 34]. 

From the above discussion, we obtain the necessary condition for the sequences uni- 
formly distributed. Moreover, if the set lim7v_;. +00 Sn given by Eq. (86) is uniformly 
distributed, Eq. (92) holds. Indeed, the Weyl criterion given by Eq. (92) is a necessary 
and sufficient condition for limjv^+oo Sn to be uniformly distributed. We give a rigorous 
proof of this fact in Appendix C. 

Let us consider a concrete example of Sn as follows: 

S N = {nAt (mod 2tt) : n G {0, 1, 2, TV}}, (93) 

where At is a positive constant. Substitution of Eq. (93) into the Weyl criterion given by 
Eq. (92) yields 

1 N I I _ im(N+l)At 
li m Y^mnAt = Um . (94) 

N^+oo N + 1 N^+oo N + 1 1 - e rmAt 

If mAt mod 2ir ^ Vm e {±1, ±2, ...}, the right-hand side of Eq. (94) is equal to zero. 
Hence, at least, letting At be written as 

At = a-, (95) 

where arbitrary positive integers p and q are coprime, and a is unity or an arbitrary pos- 
itive irrational but not transcendental number, we consider limjv^.+oo $n to be uniformly 
distributed for [0, 2ir). To derive this result, we make use of facts that tt is a transcendental 
number and a multiple of 2n is never equal to 

nAt = ncr^, (96) 
28 



where n is an arbitrary positive integer. In fact, not only At given by Eq. (95) but also a 
large variety of other time displacements lets limjv^+oc, Sn given by Eq. (93) be uniformly 
distributed for [0, 2n). However, throughout this paper, we concentrate on At given by 
Eq. (95). 

From the above considerations, we can conclude that the sequence liniAr^oo^Ar such 
that Eq. (93) is uniformly distributed modulo 2it, where At is given by Eq. (95). Moreover, 
we pay attention to the following fact. If we assume 

7i < At < 2vr, (97) 

we can consider liniA^ +00 Sn to be a pseudorandom sequence uniformly distributed for 
[0,27r). This is because a relation 2At > 2ir generates an effect of the linear congruence 
method on the sequence nAt modulo 2n for n — 0, 1, 2, .... 

We can confirm this effect by giving a concrete counter-example for < At <C ix. For 
example, taking At = 0.05, we obtain 

125At = 6.25, 126At = 6.3, 125At < 2tt < 126At, (98) 

so that (0, At,2At, 125At) form a finite arithmetic progression. Elements of the pro- 
gression do not suffer from an effect of a constant modulus 2n. When we reach at 126At, 
we observe an effect of the linear congruence method first. In this case, evidently, we 
cannot think limjy^+oo to be a random sequence. 

Using the results obtained above, we discuss an invariance of S z (t)[= L 4 (i)] under a 
change of the scale At. Here, we introduce a variable b = e _/3 . Because of < (5 < +oo, 
we obtain < b < 1. Thus, we can rewrite L 4 (t) given by Eq. (27) as 

Ut) = - 1 -(l-b)+ b -(l-- b ff(t), (99) 

oo 

f(t) = ^6 n cos(2v^t). (100) 

n=l 

In the following paragraphs, we show that the function f(t) given by Eq. (100) has a 
special scaling property. 

First, we consider a set S of pseudorandom numbers (or a sequence) distributed uni- 
formly for [0, 2n) as follows: 

S = {t n (mod 2?r) : n e {0,1,2,...}}. (101) 

Substitution of elements of S into f(t) defined in Eq. (100) yields a set of J 7 as 

?={f(t ),f(t 1 ),f(t 2 ),...}. (102) 

Here, let us consider an arbitrary real but not transcendental number s(> 1) and generate 
the following set, 

S' = {st n (mod 2?r) : n e {0,1,2,...}}. (103) 
Moreover, we construct a new set T' from S' as follows: 

r = {f(st )J(st 1 )J(st 2 ),...}. (104) 
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Then, we cannot distinguish S from S' and J 7 from J 7 ' actually on condition that both 
sets S and S' are countably infinite and their cardinalities are equal to each other, that 
is, |«S| = |«S'|. We argue this statement in the following paragraphs. 
First, we think about At given by 

At = -, - < At < it 
q 2 

where arbitrary positive integers p and q are coprime. For example, we can choose At = 
7/4 = 1.75. From the rational number At, we generate a set as follows: 

Sm,m — {mAt (mod tt) : m e {0,1,2, ...,M}}. (106) 

As explained before, if we take sufficient large M, we can consider Sm,m to be uniformly 
distributed for [0, tt). In addition, because of the effect of the constant modulus tt, we can 
regard <Sm,a< as a sequence of the pseudorandom numbers. Moreover, replacing a with 
s/n in Eq. (95) that gives the definition of At, we obtain new pseudorandom sequences 
uniformly distributed for [0,7r), 

Sm,vhm = {my/nAt (mod tt) : m e {0, 1, 2, M}} for n = 1, 2, 3, .... (107) 

Here, using the function /(t) given by Eq. (100), let us construct a set, 

J r M,At = {f(mAt (modTr)) :me {0,1,2, .,M}}. (108) 

In Eq. (100), the nth term of the series f(t) is given by b n cos (2 \/nt). Thus, a set of 
arguments \fnt of this cosine function is equivalent to S M ^ At defined in Eq. (107). 

Furthermore, taking an arbitrary real but not transcendental number s(> 1), we think 
about a change of the scale At — > sAt. This change of the scale yields the following sets: 

S MysVEAt = {srrn/KAt (mod tt) : m e {0, 1, 2, .., M}} for n = 1, 2, 3, (109) 

J~M,sAt = {f(smAt (modTr)) : me {0,1,2, ..,M}}. (110) 

In the limit of M — >■ +oo, both S M ^ At and S M ^ s ^ At become pseudorandom sequences 
uniformly distributed [0, tt), because sAt satisfies the definition of At given by Eq. (95). 
This implies that we cannot distinguish S M ^ At and S Ms ^ At actually in the limit of 
M +oo. Thus, we cannot distinguish J-M,At and J~M,sAt actually under the limit of 
M — > +oo, too. Let us describe these results as S M ^ At ~ S MjS ^ At and J~M,At — J~M,sAt- 
From these discussions, we can conclude as follows: taking the change of the scale At — > 
sAt for s > 1, we observe the scale invariance, S ~ S' and J 7 ~ J 7 ', for Eqs. (101), (102), 
(103) and (104). 

We can find similar scale invariance for S x {t)[= L^t)}. Substituting b = e~P into 
Eq. (27), we can rewrite Li(t) as 

L 1 (t) = (l-b)g(t), (111) 

oo 

g (t) = J2 bn cos(Vn + It) cos( v / ^t). (112) 

n=0 
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(105) 



Substitution of elements of S given by Eq. (101) into g(t) yields a set of 

Q = {g(t ), g(ti), g(t 2 ), ■■■}■ Then, the set Q also acquires scaling properties, which are 

similar to 7 ~ P for Eqs. (102) and (104). 

To provide pseudorandom sequences to functions L^(t) and fit) defined in Eqs. (99) 
and (100), we choose (vr/2) < At < ir. However, if we think about Li(t) and git) defined 
in Eqs. (Ill) and (112), we have to choose 7r < At < 2tt for generating pseudorandom 
sequences safely. The reason why we have to be careful in adjusting At is the difference 
of a numerical factor two between arguments of the nth terms of series f(t) and g(t), that 
is, cos(2y^nt) and cos(\/n + It) cos(-^nt). Hence, we set At = 7/2 = 3.5 for Figs. 2, 3, 4, 
5 and 6. 

In Appendix D, we examine physical transient spectra of the atom in the cavity. More- 
over, to understand the physical meanings and the scale invariance of the discrete plots 
of the trajectories of the Bloch vector, we investigate difference between continuous and 
discrete Fourier transforms for the atomic fluorescence. To confirm the scale invariance 
of the discrete plots, we examine histograms of the samples of the atomic fluorescence at 
constant time intervals. 

7 The graph of S x (t) versus the inverse of the tem- 
perature f3 for the time t such that \S z (t)\ <^ 1: the 
numerical experiments 

Figures 2, 3 and 4 suggest to us that distribution of points with \S Z \ <C 1 depends on (5 
strongly. More precisely, on the one hand, for j3 1, points of \S Z \ <C 1 are localized 
around \S X \ ~ 1. On the other hand, for j3 <C 1, points o£ ] j <C 1 are spread over 
a range of S x G [—1,1]- Figure 7 illustrates this characteristic feature of the function 
L^{t)[= S z (t)} given by Eqs. (99) and (100) well. How to draw Fig. 7 is as follows. First, 
taking < e <C 1, we look for t n (= nAt), each of which satisfies |-^4(t n )| < e. Second, we 
plot points (j3, S x (t n )) for these times t n . 

First of all, we show that the following relation holds, 

L 4 (t) = <£> cos(2 v / ^t) = 1 Vn e {1,2,3,...}. (113) 

We can derive Eq. (113) as follows. On the one hand, if we assume the right-hand 
statement of Eq. (113) holds, we can rewrite Eq. (100) as 

/(0 = f> n = A- ( 114 ) 

n=l 1 

Substitution of Eq. (114) into Eq. (99) yields L±(f) = 0. On the other hand, taking care 
of < b < 1 and -1 < cos(2 v /nt) < 1 Vn G {1,2,3,...} in Eq. (100), we can conclude 
f(t) = 6/(1 — b) if and only if the right-hand statement in Eq. (113) holds. Thus, from 
these discussions, we arrive at L 4 (t) = if and only if cos{2^/nt) — 1 Vn G {1, 2, 3, ...}. 

Therefore, we can consider the right-hand statement of Eq. (113) to be a necessary 
and sufficient condition for L 4 (t) = 0. Then, we can rewrite this necessary and sufficient 



31 



condition as 



2t = 2^71, 
2\[2t = 2n 2 n, 
2y/3t = 2n 3 7r, 

(115) 

where ni,n 2 ,n3, ... are integers. Thus, if we assume n\ ^ 0, we obtain 

V2 = ^, V3 = ^, .... (116) 

Equation (116) insists that \/2, a/3, ... are rational numbers, so that this result causes a 
contradiction. Hence, we achieve a conclusion that L±{t) = if and only if t = 0. 

However, as mentioned in Sec. 4, it is possible that |-^4(t)| <C 1 holds for some t(> 0). 
Thus, let us think about a problem whether or not we can find t, which satisfies |I>4(£)| < e 
for < 3e <C 1. Before we try a treatment of an algebraic analysis, we examine this 
problem with numerical calculations. First, we choose an arbitrary rational number as a 
time displacement n < At < 2ir. Second, taking a sequence of the time variable, t n = nAt 
for n E {0, 1, 2, N}, we construct a set, 

{L 4 (t„) :ne{0,l,2,...,7V}}. (117) 

Third, choosing < e C 1, we gather t n , each of which satisfies |-^4(t„)| < e. Because 
of Eqs. (27) and (28), we can consider L\{— S x ) to be a function of t n and j3. Fourth, 
we plot (/3,Li(t n )) for these times t n and obtain a graph of Fig. 7. The graph of Fig. 7 
represents a dependence of S x on (3 under the condition \S Z \ < e. 

When we produce the graph of Fig. 7 with numerical calculations, we have to pay 
attention to the following facts. In general, if we fix f3 and N to certain values, the 
number of t n such that |L 4 (t„)| < e decreases rapidly as e(> 0) becomes smaller. Thus, 
with fixing (3, we have to let N be larger as e(> 0) becomes smaller and approaches 
to zero. According to Eq. (117), we let L 4 (t ), L^ti), £4(^2), Li(t N ) form a set of 
samples, whose number of elements is given by (N + 1) for a certain fixed (3. Among the 
samples given by Eq. (117), the number of t n with |L4(i„)| < e decreases rapidly as (3 
becomes smaller. Thus, to let the number of t n with |L4(i n )| < e remain constant, we 
have to cause iV to be larger as f3 becomes smaller. Because of these circumstances, for 
actual numerical calculations, we introduce the following relation, 

N = N e Cl/ P, (118) 

where iVo = 619.3 and c\ = 13.37. We apply Eq. (118) to calculations with e = 7.5 x 1CT 4 
for 1.0 < (3 < 5.0. According to Eq. (118), we obtain N = 8978 for /3 = 5.0 and 
N = 396 659 904 for p = 1.0. 

However, if we apply Eq. (118) to a case of 0.5 < (3 < 1.0, iV becomes too large and 
we cannot carry out numerical calculations actually. Thus, for 0.5 < (3 < 1.0, we put 
N = Ax 10 8 as a constant and introduce an alternative relation, 

e = e e C2/ P, (119) 
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where eo = 1.875 x 10 4 and C2 = ln4. According to Eq. (119), we change e as a function 
of /3, so that we obtain e = 7.5 x 1(T 4 for (3 = 1.0 and e = 3.0 x 10" 3 for (3 = 0.5. 

Following the above prescriptions, we plot points of S x (t n )) that satisfy \S z (t n )\ < e 
for 0.5 < (3 < 5.0 in Fig. 7. For all points plotted in Fig. 7, we set At = 3.5. In Fig. 7, the 
horizontal and vertical axes are scaled logarithmically and linearly, respectively. Turning 
our eyes on Fig. 7, we observe that distinctive curves come into existence at S x — ±1 
in the limit of (3 — > +oo and their branches grow out from their trunk and spread as 
[3 becomes smaller. In Sec. 8, we analyze properties of these curves using perturbative 
techniques. 

8 The graph of S x {t) versus the inverse of the tem- 
perature /3 for the time t such that|/S z (t)| <C 1: the 
perturbative evaluation 

In this section, we examine values of the time t that satisfy |^4(t)| «C 1 in perturbation 
theory. In Sec. 4, we suggest that Fig. 7 is an appearance of incommensurate angular 
frequencies under the quasiperiodic dynamics. Then, we explain how to obtain approxi- 
mations of their ratio with continued fractions. 

In Sec. 7, we obtain a result that L A (i) = holds if and only if t — 0. Let us regard 
1/4 (i) as a power series in b from Eqs. (99) and (100). Then, we take < e « 1 and 
investigate values of the time variable t which satisfy |-^4(t)| < e. Thus, we begin the 
perturbation method with the following relation: 

Here, taking the low temperature limit j3 >> 1, that is, < 6 <C 1, we investigate Eq. (120) 
in terms of the parameter b, such that the zero, first etc., powers of b correspond to the 
zero, first, etc., orders of the perturbation calculation. 



8.1 The first and second-order perturbations 

At first, we discuss the first-order perturbation. In first order, we can rewrite Eq. (120) 

as 

b-2eb< 6cos(2t) < b. (121) 

Taking the limit of e — y +0 for Eq. (121), we obtain t — 0, 7r, 2n, 3n, .... These results are 
not important for us indeed, because they do not give us any physical meanings. 

Next, we discuss the second-order perturbation. In second order, we can rewrite 
Eq. (120) as 

< 6[1 - cos(2t)] + b 2 [l - cos(2v / 2t)] < 2e6(l + 2b). (122) 

It is very difficult for us to find all values of t, each of which satisfies Eq. (122). Thus, 
giving up our attempts to find all t with Eq. (122), we concentrate on specifying a subset 
oft such that Eq. (122). 
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Because of cos(2t) < 1 and cos(2\/2£) < 1, both the first and second terms of Eq. (122) 
have to be larger than or equal to zero. Hence, let us consider a special case where the 
following two relations hold: 

b[l - cos(2i)] = 2eb, (123) 
b 2 [l -cos(2v / 2t)] = 4e6 2 . (124) 

We have to emphasize that some t with Eq. (122) may not satisfy Eqs. (123) and (124). 
However, we neglect this possibility throughout this section, and we give all attention on 
Eqs. (123) and (124). 

Because Eq. (123) is essentially equivalent to Eq. (121), we obtain t = qn for q = 
0, 1, 2, .... Here, we remember that we obtain q = t = as a trivial root of L 4 (t) = in 
Sec. 7. Thus, we think about q = 1,2,3, ... only. Then, from Eq. (124), we obtain 

cos(2v / 2g7r) = 1 - 4e. (125) 

This implies the following fact. There are a countably infinite number of inequalities, 

(126) 
(127) 
(128) 
(129) 



\2v2~qn 


- 2pn\ 


< 


<J(4e) < 1, 


\2v2~qn 


— 4p7T 


< 


<J(4e) < 1, 


\2v2qn 


— 6p7T 


< 


<J(4e) < 1, 


\2v2~qn 


— 8j97T 


< 


<J(4e) < 1, 



where an explicit form of S(4e) is given in a next paragraph. Among the above inequalities, 
only one inequality holds for a certain positive integer p, where p and q are coprime. 

An explicit form of a function S(Ae) appearing in Eqs. (126), (127), (128) and (129) is 
given by 

<J(4e) = |arccos(l-4e)|. (130) 

Equation (130) tells us that <5(4e) is equal to a very small positive number. From Eq. (126), 
we can derive the following inequality, 

|v^-l|<^M«l. (131) 
p 2np 

Thus, we can conclude p/q ~ \/2, so that we arrive at the fact that p/q is a rational ap- 
proximate number of \/2 finally. In other words, one of the time variables with Eqs. (123) 
and (124) is given by t — qn, where p/q is a rational approximate number of \/2. 

Then, the following question arises. How small the upper bound of \y/2 — (p/q)\ is? 
How do we estimate the numerical precision of p/q as an approximate number of \/27 We 
discuss the question of accuracy later. 

From Eq. (127), we obtain 

|^-1|<^«1. (132) 
W2p Anp 
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This inequality yields a conclusion of p/q ~ 1/V2, so that we arrive at the fact that p/q 
is a rational approximate number of l/\p2. This implies that one of the time variables 
with Eqs. (123) and (124) is given by t = qn, where p/q is a rational approximate number 
of l/y/2. 

Applying a similar discussion to Eq. (128), we obtain p/q ~ y/2/3 and i = g7r. More- 
over, applying a similar discussion to Eq. (129), we obtain p/q ~ l/(2y / 2) and t = qn. 
Putting these results together, we can derive a general formula for k — 1, 2, 3, 

- ~ t = g7r for fc = 1,2,3,.... (133) 

q k 

As a result of the above discussions, we obtain the following conclusion. In the second- 
order perturbation theory, the time variable t = qn satisfies |-^4(t)| «C 1, where p/q is a 
rational approximate number of y/2/k for k = 1,2, 3, and p and g are coprime. Going 
into details, a set of the time variables t with < |Z/ 4 (t) | < e(^C 1) includes the following 
elements t = qn. First, we consider rational approximations p/q of \f2/k for k — 1, 2, 3, .... 
Second, using a very small positive value S(4e), we can give their accuracies as 

We can rewrite Eq. (134) as follows: 

|^-^|<P«1 for k = 1,2, 3, .... (135) 
k q 2knq 

In the above equation, p/q (where q > 0) represents a rational approximate number of 
y/2/k. If we put t = qn with Eq. (135), < |L 4 (f)| < e holds. Therefore, Eq. (135) 
gives an upper bound of the error resulting from this approximation. Let us estimate the 
right-hand side of Eq. (135) at 0(1/ q). Then, Eq. (135) implies that the precision of the 
rational number approximation p/q to \/2/k has to be less than 0(1/ q). 
From the above considerations, we obtain the following result: 

l^--| = — L for A; = 1,2, 3,..., (136) 
k q cq 1+u 

where p and q(> 0) are coprime and v > 1. In addition, we suppose that c is a constant 
and close to unity. Furthermore, we want to let v be as large as possible we can. This is 
because a rational approximate number p/q becomes closer to \p2/k as v(> 1) gets larger 
and larger. 

Here, we make use of the following theorem [28, 35]. For an arbitrary irrational number 
a, there exist infinite sequences p n and q n (> 0) for n > such that p n and q n are coprime 
and 

\q n a -p n \ < —. (137) 

Qn 

In other words, there exist infinitely many rational numbers p/q such that 

\a--\<\, (138) 
q q l 



35 



where p and q(> 0) are coprime. (The rigorous proof of this theorem is given in Ap- 
pendix E.) Moreover, the following fact is very useful for us in the remainder of this 
section. For an arbitrary irrational but not transcendental number a, if we suppose that 
there exist infinitely many rational numbers p/q such that 

v I 

«--<—, (139) 
1 q l qt* v ' 

the upper bound of \x is equal to two. This fact is called the theorem of Roth, which is a 
fundamental result in the Diophantine approximation [28]. Because to prove this theorem 
is beyond our purpose of this paper, we do not touch it any more. Hence, from now on, 
we consider only a rational approximate number p/q of \[2jk for k = 1, 2, 3, ... such that 

k q q 2 

and we choose t = qn for Eqs. (123) and (124). 

In the following paragraphs, we confirm the above discussions with numerical calcu- 
lations actually. By computing first some terms in a continued fraction representation of 
an irrational number, we can efficiently obtain its rational approximate number with high 
precision. Let us describe the continued fraction of an arbitrary irrational number a as 

1 

a = a -\ j 

ai + 



^| < — for k = 1,2,3,..., (140) 



1 

a 2 H 



= [a ;a 1 ,a 2 ,...], (141) 

where ao is an integer and a\, a 2 , ... are positive integers. 

The continued fractions of irrational numbers have the following properties, whose 
proofs are given in Refs. [28, 35]. 

1. If a real number a is irrational, its continued fraction expression [ao; ai, a 2 , ■■■] is 
infinite. 

2. If and only if a is an irrational solution of a quadratic equation with rational coef- 
ficients, its continued fraction expression is periodic. 

3. For an arbitrary irrational number a = [a ; a±, a 2 , a 3 , ...], let us consider infinite 
sequences of integers (po,Pi,p 2 , •••) and (?o, 9i, 92, •••) such that 

[a ;a 1 ,a 2 ,...,a n ] = —. (142) 
Obviously, from Eq. (141), we can write down p n and q n as follows: 

P-l = 1, p = Oq, p n = dnPn-l +Pn-2, (143) 

9-i = 0, g = 1, 9n = a n q n -i + 9n-2- (144) 
Then, the following relation holds, 

^<^<^<...<^<^<^. (145) 
9o 92 94 95 93 9i 
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4. For a continued fraction expression of an arbitrary irrational number 

a = [a ; a±,a 2 , a 3 , ...], let us consider a quantity obtained by including its first (n+ 1) 
terms as shown in Eq. (142). Then the following inequalities hold, 

\a-—\< — l — for n = 0,1,2,.... (146) 
q n a n+l qi 



Let us obtain a rational approximate number of y/2/k for k — 1, 2, 3, ... by truncating 
its corresponding continued fraction. We describe the obtained rational approximate 
number as p/q. Then, because of the fourth item written above, it satisfies Eq. (140). 
[Here, we pay attention to the following fact. In general, we cannot find all rational 
approximate numbers of y/2/k with Eq. (140) by using the continued fraction.] 

The continued fraction of a quadratic irrational number \[2 is given by 

v / 2=[l;2,2,2,...] = [l;2]. (147) 

In Eq. (147), a — 1 and a-i — a 2 — ... — 2 hold, so that the figure "2" appears in 
succession with period unity. From now on, we indicate the repeating block by dots as 
shown in Eqs. (147) and (149). Moreover, we introduce the following notations, 

Xv - 2 (0) = [1] = 1, 
xv - 2 (l) = [l;2] = 3/2, 
XV2 (2) = [1; 2, 2] = 7/5, 

n 

XV2 (n) = [1;0- (148) 

Because l/y/2, y/2/3 and l/(2\/2) are quadratic irrational numbers, we can write 
down them in the following continued fraction expressions, as well, 

l/y/2 = [0;1,2,2,2,...] = [0;1,2], 
V2/3 = [0; 2, 8, 4, 8, 4, 8, 4,...] = [0;2,8,4], 
1/(2^) = [0; 2, 1,4, 1,4, 1,4,...] = [0;2,i,4]. (149) 



Moreover, we obtain 



xi/A Q ) = [°1 = °» xyA 1 ) = [°; i] = 1. xi/^(2) = [o; 1, 2] = 2/3, 

n 

Xi/^W = [0;i~KX2], 

XV2/ 3 (0) = [0] = 0, X ^/s(l) = [°; 2] = 1/2, XV2/ 3 (2) = [0; 2, 8] = 8/17, 
X^2/ 3 ( n ) = [0;2,8,4, ...,a n ], 

Xi/ ( 2V2)(0) = [0] = 0, Xi/ ( 2V2)(1) = [0; 2] = 1/2, Xi/ ( 2V2)(2) = [0; 2, 1] = 1/3, 

n 

Xi/ {2 V2)H = [0;'2,l,4^...,a n ]. (150) 
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Let us consider Eq. (140) for A: = 1,2,3,4. (We neglect the cases where k — 5,6, 7, ....) 
Then, we choose rational approximate numbers of y/2, l/y/2, y/2/ 3 and 1/(2^/2) as 
follows. At first, thinking about rational approximate numbers of y/2, we dismiss x^O), 

Xv^(H) m order to obtain good accuracy. Thus, we select the following 28 rational 
numbers for the approximation of y/2, 

47 321 114 243 

X ^ (12) ~ 33 461' X ^ (13) ~ "80 782~' '"' 

, x 1 023 286 908 188 737 , . 

Y ^(39) = . (151) 

; 723 573 111 879 672 v ; 

We do not select ^^(40), y v ^(41), y v ^(42), ... because each number of digits in their 
denominators is larger than 15 and they are not tractable in the numerical calculations 
actually. In a similar way, we obtain rational approximate numbers of I /y/2, y/2/3 and 
1/(2^2) as follows: 

13 860 . . 33 461 

Xl /^ (12) ~ 19 601' Xl /V2(13) - 
_ 299 713 796 309 065 
Xl /^ (39) ~ 423 859 315 570 607' ( } 

(12) - 362 226 480 (13) - 1 492 851 361 

' ~ 768 398 401' > ~ 3 166 815 962' '"' 

58 522 759 015 841 
X ^l^ > ~ 124 145 519 261 542' ^ ' 

6930 33 461 

Xl/(2V2)( 12 ) - 1Q g 01 > Xl/(2^2)( 13 ) - g 4 g 42 ' -> 
_ fo Q \ - 2 " 713 796 309 065 

Xi/(2V2)( 39 ) " 847 7i 8 63i U1 2 U (154) 
In fact, we can verify that the following relations hold: 

|v / 2-x^(12)| ~ 3.158 x 10~ 10 < 1/(2 x 33 461) 2 , 

Iv^-X^l 3 )! ^ 5.418 x 10~ n < 1/(2 x 80 782) 2 , (155) 

and 

1 - cos(2 v / 2 x 33 461tt) ~ 2.204 x 10" 9 , 

1 - cos(2 v / 2 x 80 782tt) ~ 3.781 x 10~ 10 . (156) 

From Eqs. (151), (152), (153) and (154), we obtain 70 different positive integers as 
denominators. [Some integers appear as denominators twice or more in Eqs. (151), (152), 
(153) and (154), so that we have to avoid counting them over again.] Putting these 70 
integers and zero together, we construct a set M., whose cardinality is equal to 71, as 

M = {0,19 601,33 461,47 321,80 782,94 642,..., 
723 573 111 879 672, 

847 718 631 141 214}. (157) 
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Figure 16: A linear-log plot of (j3, S x (t)) for t = qn where q G M. given by Eq. (157) and 
2.5 < /3 < 5.0, and for t = qn where q G M given by Eq. (161) and 1.0 < < 2.5. All 
dots are drawn with changing (3 at constant interval of 0.01. We give them a diameter 
being 4/1000 of the width of the whole graph. 



Then, we obtain 71 distinct times, t = qn for q G M. According to Eqs. (27) and (28), we 
consider S x (= Li) to be a function of t and (3 and plot the points of (/3, S x (t)) for t = qn 
where q G M. and 2.5 < (3 < 5.0 in Fig. 16. Figure 16 reproduces the graph of Fig. 7 for 
2.5 < p < 5.0 well. 



8.2 The third-order perturbation 

In this subsection, we discuss the third-order perturbation. In third order, we can rewrite 
Eq. (120) as 

f(b,t) = b + b2 + b3 [b cos(2t) + b 2 cos(2 V2t) + b 3 cos(2 Vst)} . (159) 

After discussions given in Sec. 8.1, we give up our attempts to find all t with Eqs. (158) and 
(159). Thus, we concentrate on specifying a subset of t, each of which satisfies Eqs. (158) 
and (159). Obviously, the time t with Eqs. (158) and (159) has to satisfy Eq. (122), as 
well. Thus, it is possible that some of t with Eqs. (158) and (159) belong to the set of 
t = qir for q G Ai given by Eq. (157). Hence, we try finding the time t with Eqs. (158) 
and (159) from t = qn for q G M. 

In fact, we select q G M. for t = qn with Eqs. (158) and (159) on the following 
condition. Fixing the parameter for perturbation at e ~ 0.014 78 and b = e~ 2 ~ 0.1353, 
that is, (3 = 2.0, and obtaining 2efe(l + 2b + 3b 2 ) /(b + b 2 + b 3 ) ~ 0.002, we look for q such 
that 

1 - 7(e" 2 -°, qn) < 0.002 for q G M. (160) 
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As a result of numerical calculations, we find 14 integers for q satisfying Eq. (160), 



M = {0,19 601,470 832,1 607 521, 

15 994 428, 18 738 638, 768 398 401, 

10 812 186 007,21 624 372 014, 

627 013 566 048, 8 822 750 406 821, 

30 122 754 096 401,299 713 796 309 065, 

847 718 631 141 214}. (161) 

Thus, regarding S x (— L x ) as a function of t and (5 because of Eqs. (27) and (28), we plot 
the points of (/3, S x (t)) for t = qir where q G M of Eq. (161) and 1.0 < /3 < 2.5 in Fig. 16. 
We observe that Fig. 16 reproduces the graph of Fig. 7 for 1.0 < j3 < 2.5 well. 

Before closing this subsection, we make two remarks concerning our treatments of the 
Diophantine approximation. First, we indicate a property of integers q G Ai in Eq. (161). 
Remembering Eq. (140), we expect Vg G Ai to satisfy the following relation for certain 
positive integers k and q, 

\f- P -\<\- (162) 
k q q 2 

In fact, for example, for q = 19 601 G At and q = 470 832 G At, we have the following 
two relations, respectively, 

< m ,ni9 ' ( 163 ) 



'2425 19 601' 19 601 2 ' 

^ 5 '<^. (164) 



'163 101 470 832' 470 832 2 
Second, we have to point out that values of the time co-ordinate t = qn for q G Ai and 
q G At given by Eqs. (157) and (161) never fit in the distinct times t n = nAt, At = p/q 
and 7r < At < 2n discussed in Sec. 6. This implies that the graph of Fig. 7 is plotted at 
times of rational numbers and the graph of Fig. 16 is plotted at times of irrational (and 
transcendental) numbers. We cannot find an acceptable way of dealing this difference 
between Figs. 7 and 16. 



9 Discussions 

We can obtain the JCM by applying the rotating-wave approximation to a single two- 
level atom that interacts with a single mode of an optical cavity. In this paper, we 
investigate the quasiperiodicity, which the atom and the cavity field show during the time 
evolution according to the Jaynes-Cummings interaction, where the atom and the cavity 
field are initially put into a certain pure state and a mixed state in thermal equilibrium, 
respectively. 

However, for actual experiments of the cavity quantum electrodynamics (QED) in 
the laboratory, we may not observe quasiperiodicity because of lack of the rotating-wave 
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approximation. For example, it is possible that we cannot reproduce robust experimental 
results of fractal properties and scale invariance, which we discuss in this paper. 

In Ref. [36], Milonni et al. describe the interaction between a collection of two-level 
atoms and the single-mode classical electric field as the optical Bloch equations using the 
semiclassical approximation. They are given as follows: 

X(t) = -u> Y(t), 

Y(t) = u X{t) + (2/h)pE(t)Z(t), 

Z(t) = -(2/h)pE(t)Y(t), 

E{t)+u 2 E{t) = -4nNpX(t), (165) 

where (X(t),Y(t), Z(t)) represents the Bloch vector, p is the transition dipole moment, 
E(t) represents the electric field, u and u> are angular frequencies of the atoms and the 
electric field respectively, and N is a density of the atoms. 

Milonni et al. investigate the optical Bloch equations numerically and obtain the 
following results. In the rotating-wave approximation, there is no predictions of chaos. 
On the other hand, if the initial conditions let the rotating-wave approximation fail, we 
can occasionally obtain the chaotic behaviour. 

In Ref. [37], Prants et al. consider the recoil effect caused by the centre-of-mass motion 
of the atom in the cavity QED. They extend the Hamiltonian of the JCM as follows: 

H = ^ — h ^uj a a z + hujfd)a — ^-^-{a + a + a_a)) cos kfi:, (166) 

where x and p are the atomic position and momentum operators, respectively. In the 
above Hamiltonian, we consider the atom to be in a single-mode high-finesse standing- 
wave cavity, so that the dynamics of the system is sensitive to the centre-of-mass motion 
of the atom. 

Prants et al. applies the semiclassical approximation to this extended JCM. They 
assume that an expectation value of x varies in time slowly and put a certain special 
initial conditions. Then, they show that the expectation value of x obeys the following 
equation: 

x + u 2 (l — cos Vtjqr) sinrr = 0. (167) 

Moreover, they obtain the effective Hamiltonian, from which we can derive the above 
equation of motion as follows: 

1 u 2 
H = -x 2 — cj 2 cosx H [cos(x + Unt) + cos(x — fy^r)]- (168) 

In Eqs. (167) and (168), we set x = kf(x) and r = Q t. The new angular frequencies 
uj and Qn are quantities constructed from ou a , ujf and Q , and they satisfy the relation 
oj oc fi^r 1 . The Hamiltonian given by Eq. (168) represents a particle moving in the field 
of three plane waves, cos a;, cos(:r + Vt N r) and cos(a; — VL^t). It is widely known that 
this Hamiltonian induce chaotic dynamics. We emphasize that Prants et a/.'s results in 
Ref. [37] are derived under the rotating-wave approximation. In Ref. [38], Chotorlishvili 
and Toklikishvili generalize the Hamiltonian in Eq. (166) for a three level optical atom 
and discuss its chaotic dynamics. 
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In Refs. [36, 37, 38], the interactions between the atom and the cavity field are discussed 
under the semiclassical approximations. In general, it is very difficult to exactly solve the 
interaction between the two-level atom and the cavity field as a fully quantum mechanical 
model without the rotating-wave approximation. Hence, we cannot predict whether or 
not actual experiments of the cavity QED reproduce the stable results of quasiperiodicity, 
fractal properties and scale invariance. This remains to be solved in future. 

In this paper, we investigate some remarkable properties of trajectories of the Bloch 
vector, which is governed by the thermal JCM. From careful observations of the quasiperi- 
odic behaviour of the Bloch vector under the thermal JCM, we obtain novel and inter- 
esting facts, such as equivalence to the fluid dynamics, scale invariance and relation with 
the Diophantine approximation. Throughout the latter half of this paper, to examine 
quasiperiodic trajectories of the thermal Bloch vector precisely, we borrow some useful 
concepts from the number theory, which is a branch of the pure mathematics. In Sec. 6, 
we introduce a knowledge about uniform distribution of real sequences. In Sec. 8, we 
utilize methods for obtaining rational approximation of irrational numbers by means of 
their expressions of the continued fractions. 

An appearance of Shor's algorithm lets the quantum information theory attract many 
researchers' attention. We remember that Shor's algorithm is a quantum algorithm for 
solving integer factorization and discrete logarithm problems efficiently. These problems 
are regarded as important topics in the field of the number theory and the cryptography. 
Hence, progress of quantum information theory gives us new connections between ideas 
in the number theory and the quantum mechanics. 

Although the JCM was proposed and studied by researchers of quantum optics, it 
has become familiar to those of quantum information theory. The JCM has given various 
information and understandings on quantum mechanics to physicists. The authors expect 
that we can bring out a great new variety of knowledge furthermore from the JCM, 
especially, in connection with both discrete mathematics and quantum theory. 



A Physical meanings of p, p and K 

In this section, we examine the physical meanings of p(t,x), p(t,x) and K(t,x), which 
appear in Eqs. (78) and (83). We introduce the compressible fluid for explaining it to be 
equivalent to the evolution of the Bloch vector in Eq. (69). As a result of this procedure, 
we obtain these quantities, as the density of the fluid, the pressure and the external force 
per unit mass. Because the compressible fluid itself is fictitious, we need to clarify their 
physical meanings. 

First of all, using Eq. (75), we rewrite Eq. (78) as follows: 

^P + v- Vp + 0(t)p = O, (169) 

*w = v - = 2 fe| + W < 170 » 

To solve the partial differential equation (169), we separate the density of the fluid p(t, x) 
into two parts, 

p(t,x) = R(t)W(x). (171) 
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Substitution of Eq. (171) into Eq. (169) yields 



9R ^ + 4>{t) + tJ--{v ■ V)W{x) = 0. (172) 



R(t) dt ^ w W{x) 

Looking at Eq. (172), we notice the following facts. On the one hand, the first and second 
terms depend on the variable t only. On the other hand, only the third term includes 
the variable x. However, because of Eq. (74), the velocity is described as the function 
v = v(t, x), so that the third term of Eq. (172) depends on both x and t. Putting these 
considerations together and letting <p(t) be an arbitrary function, we can divide Eq. (172) 
into the following two differential equations: 

' 9 R(t) + <f>(t) = <p(t), (173) 



R(t) dt 
1 



u . {x)K v.V)W(x) = -<p(t). (174) 
Then, let us solve Eq. (174) first. Dividing W(x) into three parts as 

W(x) = X(x)Y(y)Z(z), (175) 
and taking care of Eq. (74), we rewrite Eq. (174) as follows: 

, ^ X '( X ) , s Y '(v) , s Z \*) ,s 

Vl{t ' x) J^ + V2{t ' y) Yj^ + V3{t ' z) zU^~ ip{t) ' (176) 



where 



v(t, x) 



I vi(t,x) 

v 2 (t,y) |, (177) 
V v 3 (t,z) 



and 



Vl (t,x) = j-^Z, 
V2{t > y) = L^f' 

vs(t,z) = ^[z-L 4 (t)] + L 4 (t). (178) 
Moreover, we divide Eq. (176) into the following three differential equations, 

Mt > x) W) = _v?i(t) ' (179) 

Mt,y) ^U = _v?2(t) ' (180) 

^(M)f^ = (181) 



43 



where 

<Pi(t) + <p 2 (t) + (p 3 (t) = <p(t). (182) 
Next, using Eq. (178), we rewrite Eq. (179) as 

X'(x) LAt) , , . . 

x M = - w* (t) - (183) 

The left-hand side of Eq. (183) depends only on x and the right-hand side of Eq. (183) 
depends only on t. Thus, they have to be equal to an arbitrary constant (—Ci), and we 
obtain 

xX'{x) = -dX(x), (184) 
C 1 L 1 {t) = L 1 (t)<p 1 {t). (185) 

From Eq. (184), we obtain 

X{x) = D 1 x~ Cl , (186) 
where D 1 is an arbitrary constant. In a similar way, from Eq. (180), we obtain 

Y(y)=D 2 y- c \ 
C 2 L 1 {t) = L 1 {t)ip 2 {t), (187) 

where C 2 and D 2 are arbitrary constants. 

Next, we think about Eq. (181). Using Eq. (178), we rewrite Eq. (181) as follows: 



Z(z) [L 3 (t)/L s (t)][z-L A (t)]+Lt(t)' 



The left-hand side of Eq. (188) depends only on z. Thus, we have to let the right-hand 
side of Eq. (188) not rely on the variable t, which is independent of z. This implies 
7? 3 (t) = 0. Thus, we obtain Z(z) = D 3 , where D 3 is an arbitrary constant. 
Putting these results together, we obtain 

W(x) = D 1 D 2 D 3 x- Cl y- c \ 

<p(t) = (C 1 + C 2 )^||. (189) 

Because of Eq. (171), we obtain DiD 2 D 3 ^ 0, so that we can avoid letting the density of 
the fluid p(t, x)[— R(t)W(x)] be equal to zero at every point x at every time t. Moreover, 
we have to put C\ < and C 2 < because we do not want to let W(x) be singular at x = 
and y — 0. Similarly, because we want to let W(x) converge to a finite value in the limits of 
\x\ — > +oo and \y\ — > +oo, we put C\ > and C 2 > 0. Hence, from the above discussions, 
we set C\ = C 2 = 0. Then, from Eq. (189), we obtain W(x) = W (= constant) and 
(p(t) = 0. Therefore, from Eq. (173), we arrive at 

^R(t) + (f>(t)R(t) = 0. (190) 
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Taking care of Li(0) = L 3 (0) = 1 for Eq. (27) and Eq. (170), we obtain the solution 
of Eq. (190) as 



R(t) = i?(0)exp[- f dsMs)\ 

Jo 

HtyL 3 (ty 

Therefore, we can write down p(t, x) as follows: 



(191) 



«- x)= m = tA)- (192) 

Equation (192) tells us that the density of the fluid p depends only on t, and it does not 
rely on x. 

Next, we think around the pressure p(t,x). At first, we have to carefully notice that 
we cannot determine p(t, x) only by the Navier-Stokes equations. To derive p(t, x) in the 
fluid precisely, we need an equation of state of the fluid given by classical thermodynamics. 
To let the discussion be simple, we assume the fluid to be barotropic, such as an ideal gas 
[39] . For a barotropic fluid, the pressure is a function of the density alone as 

p = f(p). (193) 

From Eqs. (192) and (193), we obtain the pressure as 

p(t,x) = f(p(t)). (194) 

Thus, we can conclude that the pressure p(t, x)[— pit)} depends only on t and it never be 
a function of x. Thus, we obtain 

Vp = 0. (195) 
Putting these discussions and Eq. (82) together, we can rewrite Eq. (83) as follows: 

dv 

— + (v-V)v = K. (196) 
at 

From Eq. (196), we can compute the external force per unit mass. Substitution of 
Eq. (74) into Eq. (196) yields an explicit form of K as 



K = 



V [L 3 (t)/L 3 (t)][z-L 4 (t)]+L 4 (t) J 



(197) 



Equations (196) and (197) mean that only the external force per unit mass K causes and 
affects the dynamics of the velocity of the fluid v(t, x). 
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B The Hamiltonian mechanics for the compressible 
fluid 



In this section, we consider how to build the Hamiltonian, which yields dynamics of the 
compressible fluid being equivalent to the Bloch vector. In Sec. 5, we show that the 
dynamics of the Bloch vector can be described with the Navier-Stokes equations for the 
inviscid compressible fluid with zero vorticity. In the following, we construct the Hamil- 
tonian that provides the Navier-Stokes equations. Moreover, we consider the reason why 
the trajectory in Fig. 1 appears to intersect itself. At the end of this section, we explain 
that the trajectory never intersects itself in the phase space. Thus, the intersections are 
the result of a projection operation and they never violate existence and uniqueness of 
solutions for the Navier-Stokes equations. 

Now, let us consider how to construct the Hamiltonian for a fluid particle [40]. First, 
we prepare a three-dimensional real vector a = (01,02,03) as a particle label at time r. 
The vector a is called a curvilinear labeling co-ordinates. They are co-ordinates fixed 
on the fluid particle. The time variable r is a proper time, which is shown by a clock 
attached to the fluid particle. Thus, the fluid particle is determined completely by (a, r). 

Let x(a, t) be the location of the fluid particle identified by a and r. Then, the vector 
a specifies the mass element of the fluid particle, so that the Jacobian determinant of a 
gives a mass-density of the fluid. The labeling co-ordinates a remain constant following 
the motion of the fluid particle. 

Here, we describe a mass element of the fluid as 

d(mass) = da 1 da 2 da 3 . (198) 

The partial differentiation of the proper time r is equal to the material derivative as 

d D d 

T = t - Tr = D-rm +v - v < (199) 

where 

„ = («,,,«,) = (* |,|). (200) 
Then, the mass-density of the fluid is given by 

9(oi,a2,a 3 ) d(a) 



d(x,y,z) d(x) 
Moreover, we define the specific volume of the fluid as 



(201) 



1 d(x) , . 

a = - = (202) 



P d(a) 

From the above notations, we obtain 

da _ d(u,y,z) | d(x,v,z) | d(x,y,w) 



dr d(a 1 ,a 2 ,a 3 ) d(a 1 ,a 2 ,a 3 ) d(a 1 ,a 2 ,a 3 ) 
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d(x) { d(u, y, z) t d(x, v, z) d(x, y, w) 



d(a) d(x,y,z) d(x,y,z) d(x,y,z) 
du dv dw 
dx dy dz 

= aV-v, (203) 



and 

»I = 4|£ = IV.., (204) 
or p p l or p 

so that we arrive at the continuity equation, 

+ P V • v = 0. (205) 

or 

Using Eq. (199), we can rewrite Eq. (205) as 

^ + (v ■ V)p + pV • v 

= f + v-M 

= 0. (206) 
Next, we define the Lagrangian as follows: 

L = jjj 'd Z a[\^T ~ E(±) - K(x,t)}, (207) 

where E{l/p) and K(x,t) represent the specific internal energy and the potential for an 
external force, respectively. We describe Hamilton's principle as 



5 J Ldr = 0, (208) 



where 5 represents an arbitrary independent variation Sx(a,r) for the fluid particle loca- 
tions. Then, we can calculate SL as follows: 

SL = -*(.,,)] 

I j III A i <? x i dE(a)Mx) dK , 

To derive the above equation, we use Eq. (202). 

Here, we introduce a convenient formula. For an arbitrary function F = F(x, t) = 
F{a, r), the following relation holds: 



d{a) 

d(5x,y,z) d(x,5y,z) d(x,y,8z) 



d 6 aF5- 

III ^ aF ^'d(a 1 ,a 2 ,a 3 ) ' d(a 1 ,a 2 ,a 3 ) ' d(ai, a 2 , a 3 ) J 
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Iff d 3 aF d ( x \ d( ' 5x,y,z " 1 + d ( x > s y> z ) + d ( x > y^ z y 

JJJ d(a) d(x,y,z) d(x,y,z) 



d{a) d(x,y,z) d(x,y,z) d(x,y,z) 

= JJJ d 3 xF{V ■ Sx) 

= - JJJ d 3 x5x ■ VF + J J F5x ■ da 

= - JJJ d 3 xSx • VF. (210) 

In the above derivation, we drop the second term as a surface integral on the two- 
dimensional boundary. 

Thus, we can rewrite Eq. (209) as follows: 



5L = - JdrJJJ d 3 xp(-^ + -Vp + VK) ■ 5x, (211) 



,d 2 x i 
p 

where we use the thermodynamical relation, 



Hence, we obtain the following equations of motion: 

8 2 t 1 

— + -Vp + VK = 0. (213) 
or 2 p 

Then, using Eqs. (199) and (200), we can rewrite the above equations as 

+ -Vp + VK = 0. (214) 
Dt p 

Consequently, we obtain the following equations of motion: 

dv 1 

— + (v ■ V)v + -Vp + VK = 0, (215) 
at p 

d -^ + V-{pv) = 0, (216) 
dE(l/p) 



on /p) +p = °- (217) 

Comparing Eq. (215) with Eq. (83), we obtain 

VK = -K, 
v = V$, 

K = -H'T-T'-lfiT-^-^ < 218 > 

and p, K, p and $ are given by Eqs. (192), (197), (194) and (80), respectively. 
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Because the Lagrangian density given by Eq. (207) is expressed as a function of (x, t), 
some might consider the number of degrees of freedom to be equal to three. However, the 
integration over the labeling co-ordinates a, that is, /// d 3 a in Eq. (207), gives an infinite 
number of degrees of freedom to the system. The Lagrangian in Eq. (207) is a sum of 
Lagrangians, each of which determines dynamics of each single fluid particle. Here, let us 
rewrite the Lagrangian with a field of the flow <p(x, t) and a field of the density p(x, t) as 



■/// 



d 3 x[\p{V V ) • (Vy?) - E(p) - K(x,t)p], (219) 



where 

d 2 E 1 dp(p) 



(220) 



dp 2 p dp 

In the derivation of the above Lagrangian, to replace labeling co-ordinates a(x, t) with 
the field of the density p(x,t) is essential. 
Then, the Hamiltonian is given by 

H = jjjd 3 x[±p(V<p)-(V<p) + E(p)+K(x,t)p]. (221) 

Hence, we can derive equations of motion as 

% = -^ = -\(V<p)-(V<p)-E'(p)-K(x,t), 
v = Vy?. (222) 

From the above results, we obtain 

-(v ■ V)w - E"(p)Vp - VK. (223) 

Because 



dv 
~dt 



E"(p)Vp = -^V P =-Vp(p), (224) 
p dp p 

we obtain the Navier-Stokes equations, 

_J = _( v .V)v--Vp + K. (225) 
at p 

In Fig. 17, the trajectory of the Bloch vector intersects itself. As explained in Sec. 4, 
for the quasiperiodic system, the trajectory of the motion never intersects itself during 
a finite time interval in the phase space because of incommensurate angular frequencies. 
The quasiperiodic trajectory becomes dense on the surface (submanifold) specified by the 
integrals of motion in the phase space. This fact seems to contradict the intersection 
of the trajectory of the Bloch vector in Fig. 17. We examine this point with numerical 
calculations. 
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Figure 17: The trajectory of the Bloch vector S(t) for < t < 5, where S(0) = (1,0,0) 
and j3 = 1.0. The trajectory intersects itself at ti = 1.644 and t 2 = 4.809. 

Figure 17 shows the trajectory of the Bloch vector S(t) for < t < 5, where S(0) = 
(1,0,0) and (3 = 1.0. Looking at Fig. 17, we notice that the trajectory intersects itself 
on (S X ,S Z ) = (0.063 72,-0.4840) at t x = 1.644 and t 2 = 4.809. However, in the phase 
space, the orbit of the motion never closes on itself at time variables ti and t 2 . We can 
confirm this fact as follows. Regarding the system as the inviscid compressible fluid with 
zero vorticity, we can understand dynamics of the system through the Hamiltonian given 
by Eq. (221). The velocity vectors (v x ,v z ) at t\ and t 2 are given numerically as 

v^tx) = -0.3632, v g (ti) = 0.2681, 

v x (h) = -0.7494, v g (t 2 ) = -0.043 48. (226) 

Hence, in the phase space, states of the system at ti and t 2 are clearly different from each 
other. 



C A proof of the Weyl criterion 

In this section, we give a rigorous proof that a real sequence is uniformly distributed 
modulo unity if and only if the Weyl criterion holds [28, 31, 32, 33, 34]. This result plays 
an important role in Sec. 6. 

Let (xo, ...,xjv) be a real sequence, which consists of (iV+ 1) real numbers. Moreover, 
for the real number Xj, let [xj] denotes the integral part of Xj, that is, the greatest integer 
not larger than Xj. Then, let {xj} = Xj — [xj] be the fractal part of Xj, or the residue of 
Xj modulo unity. Clearly, {xj} G [0, 1) holds. 

Let the counting function A(E; N+l) be defined as the number of terms Xj, < j < N, 
for which {xj} G E. Then, we define a sequence uniformly distributed modulo unity as 
follows. The real sequence (xq, ...,xn) is uniformly distributed modulo unity if for every 
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pair a, b of real numbers with < a < 6 < 1, we have 

A([a,b);N + l) , , 
lim u '- = b-a. 227 

Af^+oo N +1 

From now on, we prove that the following two statements (i) and (ii) are equivalent 
to each other. Especially, the statement (ii) is called the Weyl criterion, so that it is a 
necessary and sufficient condition of the statement (i). 

(i) A real sequence of (x ,Xi, ...,xn) is uniformly distributed modulo unity. 

(ii) For all integers / ^ 0, the following relation holds: 

1 N 

lim Y exp(2mlxi) = 0. (228) 

Before we prove (i) <^ (ii), we pay attention to the fact, exp(2irixj) = exp(2ni{xj}). 
Thus, in the remainder of this section, we consider {xj} to be equivalent to Xj without 
loss of generality. 

First, we prove (i) => (ii). At first, we suppose a real sequence (x , x±, x^) to be 
uniformly distributed modulo unity. Now, let X[a,b] be the characteristic function of the 
interval [a, b] C [0, 1), 

f \ / 1 x e[a,b], 
^ = (o xt[a,b]. (229) 

Using Eq. (229), we can rewrite Eq. (227), which defines uniform distribution of real 
sequences, as follows: 



1 N f 1 
at^+oo ivTT ? X[aM = Jo X[aM ^ dX ' (230) 



1 N /-i 

3=0 

Then, let f(x) be a step function, 

f( X ) = C X[ ,a)(x) + C a X[a,b)(x) + CbX[b,c)( X ) + - + ^X^l)^)' ( 231 ) 

which corresponds to the subdivisions, 

[0, 1) = [0, a) + [a, b) + [6, c) + ... + 1), (232) 

where cq, c a , Cb, c z are real numbers. [In Eq. (232), we divide [0, 1) into 27 parts with 
points a, b, c, z, for example. However, in fact, we can divide the interval [0, 1) into 
any number of parts.] Thus, we can consider the step function to be a superposition of 
characteristic functions with real coefficients. From Eqs. (230), (231) and (232), we obtain 



1 N r 1 
lim — — T Y,f(xj)= f(x)dx. (233) 
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Figure 18: These graphs show how to obtain the step function f(x), which approximates 
an arbitrary real function g(x) well. 

Here, let g(x) be an arbitrary real valued function, which is continuous for every value 
of x such that x G [0, 1) with g(0) = g(l). Then, 3e(> 0), we can always find a step 
function f(x) such that 

\g(x) -f(x)\ < e, 0<Vx<l. (234) 

We can understand this fact at ease by looking at Fig. 18. That is to say, f(x) is an 
approximate function of g(x), so that we can always let its accuracy be as fine as possible 
by adjusting widths of the steps and the real coefficients for the superposition. 
Then, from Eq. (234), we obtain the following inequality, 



1 JV r 1 

1 N l N ! 

< I^EJ^)-/(^]l + I^E/fe)-/ f(*)dx\ 

+ | / f(x)dx — / g(x)dx\ 
Jo Jo 

^ ^+\ J ^ Y 'tf(x j )- £ f(x)dx\. (235) 



3=0 

Thus, using Eq. (233), we can rewrite Eq. (235) as follows: 

i N 

SUD I— 



1 N r 1 
lim sup |— — Y,9{xj) ~ / g{x)dx\<2e. (236) 



3=0 

Then, remembering that e(> 0) is an arbitrary constant, we obtain 

i N 



1 N /"I 

lim — g(xj) = / q(x)dx. (237) 

-^+oo N + 1 11 Jo y ' y ' 



Let us substitute g(x) = e 2mlx for all integers I ^ into Eq. (237). Although we assume 
g(x) to be a real valued continuous function, we can put it into Eq. (237). The reason why 
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is that we can divide it into two trigonometric functions e 2mlx = cos(27r/a;) + % sin(27rZa;), 
where each of them [cos(27r/:r) and sin(27r/a;)] is real and continuous. 
Obviously, the right-hand side of Eq. (237) is equal to zero for I ^ 0, 



/' 

Jo 



e 2nilx dx 

o 



,,2nilx 



2ml 

Thus, finally, we arrive at the following relations, 

N 



0. (238) 

x=0 v ' 



lim — - — Y exp(2mlxi) = 0, (239) 

where /(^ 0) is an arbitrary integer. Hence, we conclude that (i) => (ii) is true. 

Second, we prove (ii) =>- (i). At first, we assume that the statement (ii) is true. Then, 
let us consider a function g(x), 

m 

9{x) = Yl a k exp(2TTil k x), (240) 
fc=i 



where h,h, ■■■Jm are arbitrary integers. Applying Eq. (228) to Eq. (240), we obtain 

i N 



3=0 

In the derivation of the above, we use 



1 N r 1 
lim — — T y2g(xj)= g{x)dx. (241) 
^+oo N+ 1 £^ Jo 







exp(2nil k x)dx = { ? \ k ^ °' (242) 
1 Zfc = 0. 



Next, we consider an arbitrary function f(x), which is continuous for every value 
of x such that x G [0,1) with /(0) = /(l). Then, 3e(> 0), adjusting a k and l k for 
/c = 1, 2, 3, m properly in Eq. (240), we can always let \f(x) — g(x)\ < e for < x < 1. 
This treatment means that we take a Fourier- series approximation for f(x). As a result 
of these discussions, we obtain 



1 JV r 1 

I y f(xj) — / /YxWxl 

N + ip JK 3J Jo K J 1 

+ | / g(x)dx — / /(x)dx| 
io JO 

1 N /"I 

^ 2e +I^TlE^)- J / (243) 



J"=0 

Moreover, substitution of Eq. (241) into Eq. (243) yields 

i N 

SUD I— 

N-++. 



1 * /"I 

lim sup I—— £ /(*,-) - / f(x)dx\ < 2c. (244) 
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Figure 19: These graphs explain how to choose two continuous functions fi(x) and f 2 (x) 
with Eq. (246) for the characteristic function X[a,b](x). 



Here, remembering that e is an arbitrary positive constant, we can rewrite Eq. (244) as 

N 

N^+c 



1 N r 1 
1™ T7-tE/(^)= / f(x)dx. (245) 



Now, let us take an arbitrary interval [a, b] C [0, 1). Then, 3e(> 0), we can choose two 
continuous functions fi(x) and ^(a?) such that 

/i(0) = /i(l), / 2 (0) = / 2 (l), 

fl < X[a,b] < f'2, 

I\f2{x)-h{x)]dx<~e. (246) 
Jo 

We can understand this fact at ease by looking at Fig. 19. 
On the one hand, using Eqs. (245) and (246), we obtain 

I N I N 

1 

fi(x)dx 







> / f 2 (x)dx - e 
Jo 

> J X[a,b]( x ) dx ~ e- (247) 



On the other hand, from Eqs. (245) and (246), we obtain 

I N I N 

lim sup Yr a h](xi) < lim sup f2(xj) 

1 



f 2 {x)dx 
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Jo 

< [ X[a,b]( x )dx + i- (248) 

J 



Here, remembering that e is an arbitrary positive constant and putting Eqs. (247) and 
(248) together, we arrive at 



1 N /-i 



This implies that the sequence (xq, x\, xn) is uniformly distributed modulo unity. 
Hence, we conclude that (ii) =>- (i) holds. 



D Physical transient spectra of the atom in the cav- 
ity 

In this section, we consider spectra of the atom, which develops according to the Jaynes- 
Cummings interaction with the cavity field. We discuss whether or not we can detect 
quasiperiodicity in the spectral analyses. Moreover, to understand physical meanings of 
the discrete plots of the trajectories of the Bloch vector, we investigate difference between 
continuous and discrete Fourier transforms for the atomic fluorescence. At the end of 
this section, to confirm the scale invariance for the discrete plots of the trajectories of the 
Bloch vector, we examine histograms of the samples of the atomic fluorescence taken at 
constant time intervals. 

First, we regard the time variable t as continuous one, and we define the physical 
transient spectrum of the radiation emitted by the atom as 

S(Q) = 2T [ T dh [ T dt 2 e-( r - ifi> )( r -* 1 >-( r+ifi> )( r -* 2 )(^(ti)|<T + a_|^(t 2 )), (250) 

J — CO J— CO 

where T _1 represents the filter's response time and T represents the time at which the 
measurement takes place [41, 42, 43]. The reason why the physical transient spectrum 
S(u) is given by Eq. (250) is as follows. 

Let us think about an experimental setup shown in Fig. 20 for measuring an atomic 
spectrum. We assume that we measure a state of the atom by detecting its fluorescence 
through the filter. Here, we neglect reduction of the wave function of the atom after each 
detection of its fluorescence. 

At first, let f(t) represent the elctric field of the fluorescence emitted from the atom. 
[Here, we assume that f(t) is not an operator but a function of the strength of the 
fluorescence.] Next, we suppose the filter is linear. Thus, we detect f(t) through the 
linear filter, which transforms the function f(t) into another function g(t) as 

9(t) = T{f(t)}, (251) 
where T{-} denotes a linear transformation. 
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Figure 20: An experimental setup for measuring physical transient spectra of the atom 
in the cavity. The atomic fluorescence is detected through the filter. 

Now, we think the input function f(t) to be an impulse concentrated at t — a, so that 
we describe it as f(t) = 5(t — a). Then, the impulse response of the filter is given by 

h(t;a)=T{8(t-a)}. (252) 

We can write down an arbitrary function f(t) as a superposition of the delta functions, 

/+oo 
de/(0*(*-0- (253) 
-oo 

Thus, the response of the filter for the above f(t) is given by 

9(t) = T{f(t)} 

/+oo 
dmms(t-o} 
-oo 

/+oo 
dU(0h(t;0- (254) 
-oo 

Here, we assume that the linear transformation caused by the filter is stationary, so 
that it does not depend on the time variable t. Now, let us describe the translation on the 
function of time f(t) by (— r) as f T (t) = f(t — r). Then, the fact that the transformation 
of the linear filter is stationary implies 

(T{f})r = T{f T }. (255) 

Thus, if we set h(t) = T{5(t)}, we obtain 

h(t;Z)=T{6t(t)} = ht(t) = h(t-£). (256) 

Hence, the response of the filter, whose transformation is stationary, is given by 

/+oo 
dU(0h(t-0- (257) 
-oo 
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Figure 21: The response of the filter, h(t) = 6(t)e r *. The passband width Y 1 generates 
the major contribution of the response function to the output of the filter. 



In the following, we examine properties of h(t) with giving concrete examples. 

First, we consider the case where the response of the filter is given by the Heaviside 
step function as h(t) = 9(t). Then, we obtain 

/ + OO ft 
dZf(Z)0(t-Z)= dCf(0- (258) 
-oo J — oo 

The above equation means that g(t) reflects the history of /(£) during £ e (— oo,t\. 
Because the output of the filter depends on past and current inputs but not future ones, 
we can say that the filter obtains causality. 

Second, we consider a filter whose response is given by h(t) = 9(t)e~ rt . Then, we 
obtain 

/+oo 
de/(Oe- r(t_€) 0(t-O- (259) 
-oo 

As shown in Fig. 21, setting h(t) = 9(t)e~ rt , the output of the filter depends mainly on 
the inputs /(£) for £ e (t — T -1 , t). That is to say, T" 1 corresponds to the passband width. 

Third, we think about the case where the response of the filter is given by h(t) = e~ tuJt . 
Then, we obtain 

/+0O 
dmOe~ Mt -° = F{u)e-™\ (260) 
-oo 

where F{Gi) represents the Fourier transform of f(t). The above equation implies that the 
filter extracts an angular frequency component F(u) from the input f(t) as the response. 

Putting above discussions together, we obtain a real response of an actual filter as 
follows: 

h(t,u,T) oc e^e'^^K (261) 
Thus, we can define the physical transient spectrum for the atom in the cavity as 

/+oo /'+00 
dt! / dt 2 6(T - tJe-F-^F-VeiT - t 2 ) e -(r+-)(r-t 2 ) 
-oo J ~oo 

x(i>(t!)\a + a-\4>(t 2 )), (262) 
so that Eq. (250) is derived. 
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From now on, we calculate the physical transient spectrum of the atom evolving under 
the thermal JCM according to Eq. (250). At first, we rewrite Eq. (250) for cases where 
the state of the atom is described with a density operator, 

S(u) = 2T [ T dh [ T dt2e -(T-i*)(T-ti)-(T+i*)(T-t 2 ) 

J — CO J — CO 

xIYp[C/(t 2 )p AP (0)C/t(i 1 ) (7+(7 _], (263) 

where pap(0) = Pa(0) ® pp. 

Here, as discussed in Sec. 3, we assume the initial state of the atom in the form 

Pa(0) = ^(|0) a + |1>a)(a(0|+a(1|) 

-\{\\\ <*> 

Moreover, we assume that the initial state of the cavity field obeys the Bose-Einstein 
statistic as Eq. (15), 

p p = (1 - e-P**) exp(-^wat a ). (265) 

Then, we obtain 

2 \ r] 



a-U(t 2 )p A (0) ® p?U\t l )a + = U n ° ) , (266) 



where 

V = Kofe) + M i(t 2 )]pp[4o( t i) + «oi(*i)]> (267) 

and Moo and uqi are given by Eq. (12). 

Then calculating the partial traces of the operators, we obtain 



Trp[uoo(t2)ppt*oo(*i)] = J2\. cos (\ D ( n ) + 9 2 h) + -Au 



i sm(yjD(n)+gH 2 ). 



x(l - e -P^) e -^ 



x [cos(V-D(n) + gHx) - -Aw v . 



rp r (+ N f /, \n V* 2 -/3^N Sil ^y / ^M t2 ) -^^n Sin (y^H*l) 

Tr P [Moi(^)ppMoi(ti)J = 2^#n(l-e p ) 7^= e 7^= , 

n=l ^(n) 

TrpKo(t 2 )pp4i(^i)] = 0, (268) 

where D{n) is given by Eq. (22). Here, as discussed in Sec. 3, we set Au = 0, and we 
obtain 

00 

Tr P Ko(t 2 )p pMo o(*i)] = (1 - e~^) ]T cos(|<?| yfct 2 ) cos(|(?| V^e^"- 1 ), 

n=l 

00 

Tr P K 1 (t 2 )p pMoi (t 1 )] = (l-e^)^sin(|( 7 | v ^t 2 )sin(|( 7 | v ^t 1 ) e - ta . (269) 



n=l 
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Figure 22: The spectrum S(Co) for T = 100.0, T = 0.01 and (3 = 4.0. We can observe two 
peaks at Co — 1 and Co = 

Thus, we arrive at 

Tr[f/(t 2 )p AP (0)f/t(t 1 ) (7+(J _] 

1 oo 

= o(l " e"^) Eh(lffl^) cos(b|^i)e-^^-^ 

2 n=l 

+ sin(|( ? | v ^t 2 )sin(|^| v ^t 1 )e~ /3?l H- (270) 

Here, we set h = 1. Then, replacing /3u and \g\t with /3 and t, we substitute Eq. (270) 
into Eq. (250). With tough calculations, we obtain 

S(Cb) = T [ T dt x [ T rft 2e -(r--)(T-t 1 )-(r+-)(^t 2 ) (1 _ e -/3 ) 

J —oo J —oo 

oo 

x ^2[cos(y/nt 2 ) cos( v / ^ti)e /3 + sin( v / nt 2 ) sin(v / nti)]e _/3n 



n=l 



= rE 



^ 2[ra + (r - ^) 2 ][n + (r + tCo) 2 ] 
x (e 2 ^-l)(n + r 2 + w 2 ) 

+ ( e P - l) 2 [(r 2 + w 2 - n) cos(2v^T) + 2r v / ^sin(2 v ^T)] 



(271) 



In the above equation, we assume that V and Co are in units of to and T is in units of \g\ 1 . 

Looking at Eq. (271), we notice S(—Q) = S(Co) and S't(u)) = S^u)). We show the 
spectra 5(a)) of TT = 1, T = 100.0 and /3 = 4.0, 2.0, 1.0, 0.5 in Figs. 22, 23, 24, 25, 
respectively. In Figs. 22, 23, 24 and 25, we can observe peaks of Co — y/n for n — 1, 2, 3, .... 
This is manifestation of the quasiperiodicity of the system. Moreover, if /3 becomes 
smaller, that is, if the temperature becomes higher, the number of peaks of S(Co) increases. 

Here, let us think about the sampling theorem. The sampling theorem tells us that if 
a signal has an upper angular frequency limit of Co, then we need sample points at time 
intervals being equal to or less than it /Co for enabling us to reconstruct the original signal. 
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Figure 23: The spectrum S(w) for T = 100.0, T = 0.01 and (5 = 2.0. We can observe four 
peaks at u — 1, \/2, \/3, 2. 



S 

1 r 



1 

2 



WW 



u 



CO 



Figure 24: The spectrum S(u) for T = 100.0, T = 0.01 and /3 = 1.0. We can observe 
peaks at Co — 1, \/2, \/3, 2, \/5, v^6, .... 
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Figure 25: The spectrum S(u) for T = 100.0, T = 0.01 and (3 = 0.5. We can observe 
peaks at u — 1, v^2, v^3, 2, v^5, v^6, \/7, 2v^2, 3, .... 
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Thus, to obtain information about the component of u — ^Jn, we have to take samples 
of the signal at time interval being equal to or less than At = 7r / \fn. On the other hand, 
from the discussions about the scale invariance and the pesudorandom sequence uniformly 
distributed for [0, 2n) in Sec. 6, we obtain n < At < 2tt in Eq. (97). Thus, for discritizing 
the dynamics of the Bloch vector, we have to set tt < At < 2ir. In fact, in Figs. 2, 3, 4, 
5, 6 and 7 in Sec. 3, we take At = 3.5 to convert the continuous trajectory of the Bloch 
vector to a discrete-time sequence. 

From these discussions, we understand that discretization of the dynamics of the 
Bloch vector with the time interval 7r < At < 2tt destroys information of all components 
of ui — \fn for n = 1,2,3,..., so that we cannot recover the original trajectory. This 
implies that the discretization with 7r < At < 2n erases the history of the trajectory of 
the Bloch vector completely. Because of this effect, figures consisting of the discrete-time 
sequence of the Bloch vector acquire the scale invariance and other fractal properties. 

In the following, we discuss the Fourier transform of the discrete-time sequence ob- 
tained along the trajectory of the Bloch vector with the time interval tt < At < 2tt. Let 
us define the discrete Fourier transform of the fluorescence of the atom in the cavity as 
follows: 

1 N-1N-1 

c t = ^EE e-^- l ^ N Tr[U(l 2 At) P M0)UHhAt)a + a_]. (272) 
iV i 1= o i 2 =o 

Looking at the above equation, we notice Ck+N = so that we concentrate only on 
{cfc : k = 0, 1, N — 1}. With tough calculations, we obtain 



Cfc 



oo i a 

V - e~ nf} sinhf-1 

^iV 2 [cos(2A;7r/iV) - cos^At)] 2 V 



x 



8 2krr 
cosh(|)[(l - cos(-^) cos(VriAt))(l - cos(2A;7r) cos(y/nNAt)) 



2kir 

- sin(2A;7r) sin(-^) sin(y/nAt) sm(y/nNAt)] 
8 2kir 

+ sinh(|) [cos(^p) - cos( v / nAt)] cos(y/n(N - 1) At) 
x [cos(2A;7r) - cos(^NAt)] 



(273) 



In Figs. 26 and 27, we plot {\c k \ : k = 0,1,2,...,N - 1} for 8 = 1.0, At = 3.5, 
N = 4000 and 8 = 1.0, At = 4.5, TV = 4000, respectively. Comparing Figs. 26 and 27, 
the Fourier transforms of At = 3.5 and At = 4.5 are obviously different from each other. 
This observation seems to insist that the scale invariance does not hold. That is, under 
a rescaling from At = 3.5 to At = 4.5, the figure consisting of the discrete-time sequence 
of the Bloch vector changes. However, this discussion is not true. 

In fact, it is not significant to apply the Fourier analysis to behaviour of the discrete- 
time sequence of the Bloch vector. This is because we take a time interval tt < At < 
2ir and we cannot recover entire information about components of Co — \fn for n = 
1,2, ...,N — 1. That is to say, a set of samples taken at time interval At lose the past 
history of the trajectory. When we apply the discrete Fourier transform to iV samples of 
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Figure 26: The discrete Fourier transform of the Bloch vector for /3 — 1.0, At — 3.5 and 
N = 4000. 
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Figure 27: The discrete Fourier transform of the Bloch vector for /3 — 1.0, At — 4.5 and 
N = 4000. 
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Figure 28: The histograms of sets given by Eq. (275) for At = 3.5, f3 = 1.0, N = 500 000 
and At = 4.5, f3 = 1.0, N = 500 000 with a thin solid curve and thick dashed curve, 
respectively. In the graphs, we let the width of bins in the horizontal axis be equal to 
0.0025. 

a signal f(t), 

/(0), /(At), /(2At), f((N - l)At), (274) 

we assume that they are samples arranged in chronological order. However, the samples 
of the trajectory of the Bloch vector lose their past history completely, so that N samples 
cannot form a sequence of the chronological order shown in Eq. (274). 

We can only regard the samples shown in Eq. (274) as a set 
{/(nAt) : n = 0, 1,2, ...,N — 1}, which lose the information about their order. Here, 
thinking about the scale invariance At — > sAt for s > 1, we want to show that we cannot 
distinguish two sets of samples, {/(nAt) : n = 0,1,2, N — 1} and {f(nsAt) : n = 
0, 1, 2, N - 1}, under the limit of iV -> oo. 

Thus, we analyze the above two sets of the samples in the following way. First, let us 
consider a histogram of the following set, 

{Tr[U(lAt)p AP (0)U j (lAt)a + a_} : I = 0, 1, 2, N - 1}. (275) 

Second, we prepare two sets of Eq. (275) for two different time interval, At and sAt. If 
we cannot distinguish their histograms under the limit of N — > oo, we can conclude that 
the scale invariance holds under At — > sAt. 

In Fig. 28, we plot histograms of sets given by Eq. (275) for At = 3.5, (5 = 1.0, 
N = 500 000 and At = 4.5, (3 = 1.0, N = 500 000 with a thin solid curve and a thick 
dashed curve, respectively. In Fig. 28, we let the width of bins in the horizontal axis be 
equal to 0.0025. Looking at Fig. 28, we can conclude that histograms of At = 3.5 and 
At = 4.5 are quite similar and we can hardly distinguish them. Thus, we can confirm the 
scale invariance from Fig. 28. 
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E Rational approximations of irrational numbers 
and their precisions 

In this section, we prove two theorems, which are related to the Diophantine approx- 
imation [28, 35]. From these theorems, we can estimate precisions given by rational 
approximations of irrational numbers, as utilized in Sec. 8. 

Theorem 1. (Dirichlet's approximation theorem) 

Let a be an arbitrary positive number, and let Q be an arbitrary integer with Q > 1. 
Then, there is a rational number p/q such that 

I«--I<4y ( 276 ) 
where p and q are coprime and < q < Q. 

Proof. 

We prove the above theorem by the pigeon-hole principle. On the one hand, using a, 
we put (Q + 1) points (pigeons) in the interval [0, 1] as follows: 

0, a - [a], 2a - [2a], (Q - l)a - [(Q - l)a], 1. (277) 

On the other hand, we split the interval [0, 1] into Q subintervals (pigeon-holes), each of 
length l/Q, as follows: 

[0, l/Q), [l/Q, 2/Q), [(Q - 2)/Q, (Q - 1)/Q), [(Q - 1)/Q, 1]. (278) 

Then, two or more pigeons have to belong to a particular pigeon-hole. Thus, one of 
the following two statements has to hold true. The first one is that there exist integers j 
and k such that 

\(ka - [ka]) - (ja - \ja})\ < l/Q, (279) 
where < j < k < Q. The second one is that there exists an integer k such that 

\(ka- [ka]) - 1| < l/Q, (280) 

where < k < Q. 

In the case where Eq. (279) holds, we define integers p and q as 

q = k-j, p=[ka]-[ja]. (281) 

Then, from Eq. (279), we obtain 

\qa-p\<l/Q. (282) 

This implies that we achieve Eq. (276). In the case where Eq. (280) holds, we define 
integers p and q as follows: 

q = k, p=[ka] + l. (283) 
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Then, from Eq. (280), we obtain 



\qa-p\<l/Q. (284) 

This implies that we arrive at Eq. (276). 

Thus, we conclude that we have proved the theorem. 

Theorem 2. 

Let a(> 0) be an arbitrary irrational number. Then, there are two infinite integer 
sequences, p n and q n (where n > 0), such that 

1 , X 

\q n a -p n \ < — , (285) 
q n 

where p n and q n are coprime and q n > 0. In other words, there are infinitely many rational 
numbers p/q, where p and q are coprime and q > 0, such that 

\a-^\<\. (286) 

q q 2 



Proof. 

We prove this theorem by assuming the opposite. At first, let us suppose that there 
are only a finite number of p n /q n , say Pi/q±, £>7v/<?7V with Eq. (285). Then, we 

consider the following real number p, 

P = min -Pn|- (287) 

l<n<N 

Because of Theorem 1, there exists p/q such that 

|a--|<4y 0<g<Q, (288) 

q qQ 

where Q is an integer and larger than unity. Moreover, we can choose an arbitrary integer 
greater than unity as Q. Thus, we set Q > 1/ p. Then, the following relation holds, 

k«- P |<i i<P- (289) 

This result contradicts Eq. (287), which defines p as the minimum number. Thus, using 
proof by contradiction, we can conclude Theorem 2. 
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